[FreeGeoCZ] RE: pajtn - gdalnumeric.LoadFile
MICHÁLEK Jan Mgr. <[email protected]>
| Newsgroups | gmane.comp.gis.grass.cz |
|---|---|
| Message-ID | <[email protected]> |
Tak už vim čim to je, ten muj tiff je jen jedno pásmo. Co s tim ale... From: [email protected] [mailto:[email protected]] On Behalf Of MICHÁLEK Jan Mgr. Sent: Tuesday, May 07, 2013 1:06 PM To: Svobodná geoinformační infrastruktura ([email protected]) Subject: [FreeGeoCZ] pajtn - gdalnumeric.LoadFile Ahojte pánové, mám tady, asi triviální prblm, se kterym ale neumim pohnout. Potřeboval jsem ořezat větší množství rastrů v tifu podle brambory v šejpu. Vygůglil jsem tady todle http://geospatialpython.com/2011/02/clip-raster-using-shapefile.html Na vzorových datech mi to v osgeo4w běhá uspokojivě. Data, který tim ale chci procedit to neskousne. Problém vzniká tady ------------------------ # Load the source data as a gdalnumeric array srcArray = gdalnumeric.LoadFile(raster) # Also load as a gdal image to get geotransform (world file) info srcImage = gdal.Open(raster) geoTrans = srcImage.GetGeoTransform() # Create an OGR layer from a Field boundary shapefile field = ogr.Open("%s.shp" % shp) lyr = field.GetLayer(shp) poly = lyr.GetNextFeature() # Convert the layer extent to image pixel coordinates minX, maxX, minY, maxY = lyr.GetExtent() ulX, ulY = world2Pixel(geoTrans, minX, maxY) lrX, lrY = world2Pixel(geoTrans, maxX, minY) # Calculate the pixel size of the new image pxWidth = int(lrX - ulX) pxHeight = int(lrY - ulY) clip = srcArray[:, ulY:lrY, ulX:lrX] Na tom clipu se mě to zrůcá, protože pole z původních dat je trojrozměrný a pole z mých dat jen dvou. Když tam tu první dvojtečku zruším, tak to doběhne, ale výsledek je podivně pomršen. Nepovedlo se mě dohledat, co mi tam dělá tu paseku. Počítám, že tady bude někdo, kdo tomu rozumí. Nemohl by mě ten někdo navýst na řešení. Díkec Jelen _______________________________________________ FreeGeoCZ mailing list [email protected] http://mailman.fsv.cvut.cz/mailman/listinfo/freegeocz