Le samedi 28 février 2015 18:29:01, Yuta Sato a écrit : > Dear Even Rouault: > > Thanks for you quick response. > Unfortunately, I do not know how to pipe using python script. > I just know to use as follows: > > for x, y in zip(lons, lats): > vals = subprocess.check_output([gdallocationinfo, '-valonly', > '-l_srs', 'wgs84', image, str(x), str(y)]) > > The input lons and lats are in the text file.
As I suggested, you'd better extract the gist of http://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/val_at_coord.py So something like from osgeo import gdal from osgeo import osr ds = gdal.Open( image ) srs = osr.SpatialReference() srs.ImportFromWkt(ds.GetProjection()) srs_wgs84 = osr.SpatialReference() srs_wgs84.SetFromUserInput('WGS84') ct = osr.CoordinateTransformation(srs_wgs84, srs) for lon, lat in zip(lons, lats): (X, Y, _) = ct.TransformPoint(lon, lat) geomatrix = ds.GetGeoTransform() (success, inv_geometrix) = gdal.InvGeoTransform(geomatrix) x = int(inv_geometrix[0] + inv_geometrix[1] * X + inv_geometrix[2] * Y) y = int(inv_geometrix[3] + inv_geometrix[4] * X + inv_geometrix[5] * Y) vals = ds.ReadAsArray(x,y,1,1) > > Can you help me? > <https://plus.google.com/u/0/106381021255842353792?prsrc=4> > > On Sun, Mar 1, 2015 at 2:05 AM, Even Rouault <[email protected]> > > wrote: > > Le samedi 28 février 2015 17:51:47, Yuta Sato a écrit : > > > I have a list of longitudes and latitudes for which I extract the pixel > > > values using gdallocationinfo.exe one by one for each pairs of > > > longitudes and latitudes using for loop. > > > > > > How I simultaneously get pixel values for all the list of longitudes > > > and latitudes? > > > My list is large, so I want faster retrieval. > > > I am using python and outputting the pixel values into text file. > > > Is it possible? > > > > Yuka, > > > > you can pipe into gdallocationinfo several tuples of coordinates: > > > > $ printf "0 0\n0 1\n" | gdallocationinfo byte.tif > > > > Report: > > Location: (0P,0L) > > > > Band 1: > > Value: 107 > > > > Report: > > Location: (0P,1L) > > > > Band 1: > > Value: 115 > > > > But I think it would be easier/faster to just do that fully in Python. > > The following script has likely everything you need to do that: > > > > http://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/val_at_coord.py > > > > Even > > > > > > -- > > Spatialys - Geospatial professional services > > http://www.spatialys.com -- Spatialys - Geospatial professional services http://www.spatialys.com _______________________________________________ gdal-dev mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/gdal-dev
