Tengo una Línea, definida como 2 puntos en un sistema de coordenadas cartesianas.
Shapely permite proyectar esta línea con una transformación afín si se especifica una transformación.
from shapely.geometry import LineString
from shapely.affinity import affine_transform
coef_6_transform = [0.10000036968576852, 0.0, 0.0, -0.1000006060605726, 350482.8259, 5617200.9533]
# obtained from a .geotif raster
def to_Shapley_Line(line):
p0, p1 = line
A = Point(p0[0],p0[1])
B = Point(p1[0],p1[1])
return LineString([A, B])
def vectorize_line(line):
return affine_transform(to_Shapley_Line(line), coef_6_transform)Funciona perfectamente pero no se ajusta a mis necesidades. Exportar un archivo como .geotiff y recargarlo es demasiado lento.
Por lo que yo entiendo, es el coef_6_transform el mapeo de a desde un mapa lineal (también conocido como sistema de coordenadas catesianas) al SIR definido.
Si proporciono el CRS y una ubicación, ¿puedo obtener esta transformación sin extraerla de un .geotiff ?
O
¿Cuál es el nombre de mi matriz local, que empieza en 0? Me gustaría et con geopandas .set_crs() después pude reproyectar con .to_crs
ejemplo para aclarar :
from shapely.geometry import LineString, Point, Polygon
from shapely.affinity import affine_transform
from owslib.wms import WebMapService
import rasterio
import geopandas as gpd
from matplotlib import pyplot as plt
from rasterio.plot import show
# Geotiff Extraction; i want to extract as a jpeg instead, to save download and computation Time
WMS_url_hist = "https://www.wms.nrw.de/geobasis/wms_nw_hist_dop?"
wms_hist = WebMapService(WMS_url_hist)
crs = "EPSG:25832"
shp_geom = Polygon(((351511.3832, 5616978.270700001), (351552.41784617, 5616997.140201801),
(351636.0678615501, 5617035.606090089), (351800.7549000001, 5617111.3363),
(351804.7563000001, 5617104.792899999), (351837.5250559599, 5617051.2122095),
(351888.6184044101, 5616967.669299701), (351621.72937412, 5616845.93172588),
(351619.2559, 5616849.2452), (351605.7223457001, 5616865.068229449),
(351586.34850795, 5616887.719544729), (351568.5424999999, 5616908.5378),
(351554.46985843, 5616925.705787269), (351543.1406684901, 5616939.5276307),
(351529.7399, 5616955.877399999), (351511.3832, 5616978.270700001)))
#print("shp_geom ",shp_geom )
bounds = shp_geom.bounds
#results in :
#bb order (minx, miny, maxx, maxy)
bounds = (351511.38320000004, 5616845.93172588, 351888.6184044101, 5617111.3363)
img = wms_hist.getmap(layers=["nw_hist_dop_2013"],
styles=["default"],
srs=crs,
bbox=bounds,
#size=(3774,2654),
size=(3000, 2500), # The original resolution is 10cm
format='image/tiff',#'image/jpeg',#Here word be a change for the other Case
transparent=True
)
out = open("Test.tif", 'wb')
out.write(img.read())
out.close()
#This line is extracted from the picture above
line = ((1068, 1908), (1351, 1757))
def to_Shapley_Line(line):
p0, p1 = line
A = Point(p0[0],p0[1])
B = Point(p1[0],p1[1])
return LineString([A, B])
# i would like to infer the right 6 parameter Transform (coeff_6_transform) from the length of the bounding Box
src = rasterio.open("Test.tif")#.jpeg
transform = src.transform
coef_6_transform = [element for tupl in transform.column_vectors for element in tupl]
print("coef_6_transform",coef_6_transform)
# results in :
# coef_6_transform = [0.1257450681366632, 0.0, 0.0, -0.10616182964801991, 351511.38320000004, 5617111.3363]
# the 2nd and 3rd parameter need to be 0, as i observed
# the 5th and 6th parameter are the minx and maxy of the bounding box
# the 1st and 4th parameter are due to the extraction of a 10 cm spatial resolution image wich a size of 3000,2500
# (3774,2654) would be more accurate and resilt in theese coefficients beeing closer to 10, but not exactly 10, unfortunatly
# is there a way of calculating the 1st and 4th paramter, just from the bounding Box and the 10cm resolution ?
def vectorize_line(line):
#transforms the line into the desired crs
return affine_transform(to_Shapley_Line(line), coef_6_transform)
transformed_line = vectorize_line(line)
fig, ax = plt.subplots()
show(src, ax=ax)
gpd.GeoSeries(transformed_line).plot(ax=ax, color = "red", linewidth = 1)
plt.savefig("Test.png")
plt.close()
