Date   

Re: rasterio opens file from AWS S3 bucket on local machine, but can't find file when deployed to Google App Engine

Judson Buescher
 

Hi Sean,

The only way that I’m installing the package is through my requirements.txt which I include in my deploy folder. I have

 

rasterio==1.1.5

 

in the requirements.txt. For the credentials I’ve tried a couple different ways. I’ve tried setting them using

 

os.environ['GS_SECRET_ACCESS_KEY'] = secret_access_key

os.environ['GS_ACCESS_KEY_ID'] = access_key

 

I’ve also tried

 

session = boto3.Session(aws_access_key_id=access_key,

                        aws_secret_access_key=secret_access_key)

 

and then using that session object when I open the raster file. Was there any other info you’d like?

Thanks a ton for the help,
Judson

 

From: <main@rasterio.groups.io> on behalf of "Sean Gillies via groups.io" <sean.gillies@...>
Reply-To: "main@rasterio.groups.io" <main@rasterio.groups.io>
Date: Saturday, August 8, 2020 at 10:43 PM
To: "main@rasterio.groups.io" <main@rasterio.groups.io>
Subject: Re: [rasterio] rasterio opens file from AWS S3 bucket on local machine, but can't find file when deployed to Google App Engine

 

Apologies for the late reply. There may be a bug, yes, but there are more likely sources of trouble. Can you explain more about how you are installing the rasterio package on GCP and how you are configuring your AWS credentials?

 

On Thu, Jul 23, 2020 at 7:06 AM <judson.buescher@...> wrote:

Hello,

 

I'm trying to run a python flask app, that opens a file hosted on an S3 bucket and analyzes the data. When I deploy the flask app on my local host, it runs like expected. Once I deploy the app to the cloud, it suddenly cannot find the file. I feel like this might be a bug, but i'm too afraid to post in the github. I'll attach the relevant snippet of code, and the error I'm receiving in the GCP Logs. I apologize for terrible formatting as I'm unaware of how to do it properly in this forum.
Here's the code: 

session = boto3.Session(aws_access_key_id=access_key, aws_secret_access_key=secret_access_key)    

def start_process(map, lon, lat):    
    with rasterio.Env(AWSSession(session)):        
        return getCoordinatePixel(map,lon,lat)  url = 's3://metcon-test-bucket/high_plus.tif'  

app = Flask(__name__) @app.route('/api')
def api_id():    
    lat = int(request.args['lat'])    
    lon = int(request.args['lon'])          
    return jsonify(start_process(url, lon, lat))  

if __name__ == '__main__':    
    app.run()

 and here's the error:



--

Sean Gillies


Re: window_transform doesn't produce accurate coordinate transforms

pranay1117@...
 

Hi Sean, thanks for your reply - that was indeed the source of the issue on my end. I was able to fix it!

--
Best, 
Pranay


Re: Reprojecting output image to input image crs where the image sizes differ

Sean Gillies
 

Hi,

On Tue, Jul 21, 2020 at 11:52 PM <ashnair0007@...> wrote:
I have a geotiff (size= 4490 x 7341) that I cut into 256 x 256 crops to do an operation. Once the operation is done, I merge them back together. Since the original dimensions aren't divisible by 256 the resulting image will be of size 4608 x 7424. However the output tif and the input tif are of the same location, so they need to be geo-referenced the same. Is it possible to do so?

If your original raster's georeferencing did not involve rotation, and if the merged raster has the same origin (upper left corner) as the original, you can reuse the original raster's transform property. If the merged raster has a different origin, you must compute a new transform using https://rasterio.readthedocs.io/en/latest/api/rasterio.transform.html#rasterio.transform.from_bounds or https://rasterio.readthedocs.io/en/latest/api/rasterio.transform.html#rasterio.transform.from_origin.

Sorry about the late reply. I hope it is still useful.

--
Sean Gillies


Re: window_transform doesn't produce accurate coordinate transforms

Sean Gillies
 

Hi,

On Thu, Jul 23, 2020 at 7:06 AM pranay1117 via groups.io <pranay1117=yahoo.co.in@groups.io> wrote:
Hi, I'm new to rasterio and working with geographical data. So far, it's been a great experience working with this library - but I've run into an issue when working with the Window module. My use case is as follows: 
- I obtain 2 bounding coordinates from the user (top-left and bottom-right) and convert them to (row, col) using `source.index(lon, lat)`.
- I create a window using these bounds and read the dataset from that. Also, I get the appropriate transform for this window.
selected_window = Window.from_slices((row0, row1), (col0, col1))
data = src.read(1, window=selected_window)
affine_transform = src.window_transform(selected_window)

Afterwards, I am using this `data` numpy array to select various points in these bounds and plot them on a map. However, the (lat, long) coordinates obtained using the affine_transform above are not very accurate. When I run the following:
print("Lat/Lon bounds: {}".format(src.window_bounds(selected_window)))
print("Transform at (0,0): {}".format(affine_transform * (0,0)))
r, c = data.shape
print("Transform at {}: {}".format((r-1, c-1), affine_transform * (r-1, c-1)))

The output is as follows:
Source lat/lon bounds: (72.45 20.69, 82.39, 28.78)
Transform at (0,0): (72.42, 28.78)
Transform at (970, 1193): (80.53, 18.841)

As is evident from the output, the latitude must lie between (20.69,28.78) but the transform value 18.841 at (970, 1193) is not in that interval. 

I saw some earlier posts in this group that suggested updating to version 1.1.5. I've done that but the issue persists. 
Any help will be greatly appreciated. Thanks!

Your usage of the affine transform

    affine_transform * (r-1, c-1)

is incorrect. You must swap the order of the tuple indexes to get the correct results because mathematically speaking

    (x', y') = affine_transform * (x, y)

and row is y and column is x. Rasterio has a method that will help keep this straight: https://rasterio.readthedocs.io/en/latest/api/rasterio.transform.html#rasterio.transform.xy.

The ordering of items in (row, col) versus (x, y) is a constant source of confusion. I hope this reply hasn't come too late and that it is useful.

--
Sean Gillies


Re: rasterio opens file from AWS S3 bucket on local machine, but can't find file when deployed to Google App Engine

Sean Gillies
 

Apologies for the late reply. There may be a bug, yes, but there are more likely sources of trouble. Can you explain more about how you are installing the rasterio package on GCP and how you are configuring your AWS credentials?

On Thu, Jul 23, 2020 at 7:06 AM <judson.buescher@...> wrote:

Hello,

 

I'm trying to run a python flask app, that opens a file hosted on an S3 bucket and analyzes the data. When I deploy the flask app on my local host, it runs like expected. Once I deploy the app to the cloud, it suddenly cannot find the file. I feel like this might be a bug, but i'm too afraid to post in the github. I'll attach the relevant snippet of code, and the error I'm receiving in the GCP Logs. I apologize for terrible formatting as I'm unaware of how to do it properly in this forum.
Here's the code: 

session = boto3.Session(aws_access_key_id=access_key, aws_secret_access_key=secret_access_key)    

def start_process(map, lon, lat):    
    with rasterio.Env(AWSSession(session)):        
        return getCoordinatePixel(map,lon,lat)  url = 's3://metcon-test-bucket/high_plus.tif'  

app = Flask(__name__) @app.route('/api')
def api_id():    
    lat = int(request.args['lat'])    
    lon = int(request.args['lon'])          
    return jsonify(start_process(url, lon, lat))  

if __name__ == '__main__':    
    app.run()

 and here's the error:




--
Sean Gillies


Re: returning subset of raster data via API

Sean Gillies
 

On Wed, Aug 5, 2020 at 5:41 AM Tom Kralidis <tom.kralidis@...> wrote:

Per https://gdal.org/user/virtual_file_systems.html#drivers-supporting-virtual-file-systems and confirmed by Even R, GDAL's NetCDF driver cannot write to /vsi*. Perhaps a workaround can be to drop to a tempfile approach when writing NetCDF.

Thanks again for the clarifications.


You're welcome. It's unfortunate that we don't have a more informative message from GDAL that we can raise to users. "No such file or directory" is too opaque.

--
Sean Gillies


Re: returning subset of raster data via API

 

Per https://gdal.org/user/virtual_file_systems.html#drivers-supporting-virtual-file-systems and confirmed by Even R, GDAL's NetCDF driver cannot write to /vsi*. Perhaps a workaround can be to drop to a tempfile approach when writing NetCDF.

Thanks again for the clarifications.


Re: returning subset of raster data via API

Sean Gillies
 

Hi Tom,

On Mon, Aug 3, 2020 at 5:32 PM Tom Kralidis <tom.kralidis@...> wrote:

(sorry, just getting back to this!)

Thanks Sean. Indeed the MemoryFile strategy works as expected:

import sys

import rasterio
from rasterio.io import MemoryFile

with rasterio.open(sys.argv[1]) as src:
    with MemoryFile() as memfile:
        profile = src.profile
        with memfile.open(**profile) as dest:
            dest.write(src.read())

Are there any known issues with reading/subsetting, and writing back a NetCDF via MemoryFile? Traceback:

(pygeoapi) (base) tomkralidis@aruba:~/Dev/pygeoapi/pygeoapi-tomkralidis$ python dd.py tests/data/CMIP5_rcp8.5_annual_abs_latlon1x1_PCP_pctl25_P1Y.nc 
Traceback (most recent call last):
  File "rasterio/_io.pyx", line 1135, in rasterio._io.DatasetWriterBase.__init__
  File "rasterio/_err.pyx", line 205, in rasterio._err.exc_wrap_pointer
rasterio._err.CPLE_OpenFailedError: Unable to create netCDF file /vsimem/fa663a74-d548-494d-a2ab-880bdff831c4/fa663a74-d548-494d-a2ab-880bdff831c4. (Error code 2): No such file or directory .

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "dd.py", line 9, in <module>
    with memfile.open(**profile) as dest:
  File "/Users/tomkralidis/opt/miniconda3/lib/python3.7/site-packages/rasterio/env.py", line 382, in wrapper
    return f(*args, **kwds)
  File "/Users/tomkralidis/opt/miniconda3/lib/python3.7/site-packages/rasterio/io.py", line 136, in open
    nodata=nodata, sharing=sharing, **kwargs)
  File "rasterio/_io.pyx", line 1139, in rasterio._io.DatasetWriterBase.__init__
rasterio.errors.RasterioIOError: Unable to create netCDF file /vsimem/fa663a74-d548-494d-a2ab-880bdff831c4/fa663a74-d548-494d-a2ab-880bdff831c4. (Error code 2): No such file or directory .

Versions:

>>> import rasterio
>>> rasterio.__version__
'1.1.5'
>>> rasterio.gdal_version()
'3.0.4'

Thanks

..Tom

I don't use netcdf files with GDAL often and have never tried to write a netcdf file to vsimem. It could be that there is a combination of creation or config options required that I don't know about. I see some mention about platform requirements in https://github.com/OSGeo/gdal/pull/786 but am not sure what to make of that.

--
Sean Gillies


Re: returning subset of raster data via API

 

(sorry, just getting back to this!)

Thanks Sean. Indeed the MemoryFile strategy works as expected:

import sys

import rasterio
from rasterio.io import MemoryFile

with rasterio.open(sys.argv[1]) as src:
    with MemoryFile() as memfile:
        profile = src.profile
        with memfile.open(**profile) as dest:
            dest.write(src.read())

Are there any known issues with reading/subsetting, and writing back a NetCDF via MemoryFile? Traceback:

(pygeoapi) (base) tomkralidis@aruba:~/Dev/pygeoapi/pygeoapi-tomkralidis$ python dd.py tests/data/CMIP5_rcp8.5_annual_abs_latlon1x1_PCP_pctl25_P1Y.nc 
Traceback (most recent call last):
  File "rasterio/_io.pyx", line 1135, in rasterio._io.DatasetWriterBase.__init__
  File "rasterio/_err.pyx", line 205, in rasterio._err.exc_wrap_pointer
rasterio._err.CPLE_OpenFailedError: Unable to create netCDF file /vsimem/fa663a74-d548-494d-a2ab-880bdff831c4/fa663a74-d548-494d-a2ab-880bdff831c4. (Error code 2): No such file or directory .

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "dd.py", line 9, in <module>
    with memfile.open(**profile) as dest:
  File "/Users/tomkralidis/opt/miniconda3/lib/python3.7/site-packages/rasterio/env.py", line 382, in wrapper
    return f(*args, **kwds)
  File "/Users/tomkralidis/opt/miniconda3/lib/python3.7/site-packages/rasterio/io.py", line 136, in open
    nodata=nodata, sharing=sharing, **kwargs)
  File "rasterio/_io.pyx", line 1139, in rasterio._io.DatasetWriterBase.__init__
rasterio.errors.RasterioIOError: Unable to create netCDF file /vsimem/fa663a74-d548-494d-a2ab-880bdff831c4/fa663a74-d548-494d-a2ab-880bdff831c4. (Error code 2): No such file or directory .

Versions:

>>> import rasterio
>>> rasterio.__version__
'1.1.5'
>>> rasterio.gdal_version()
'3.0.4'

Thanks

..Tom


Re: How to extract DATA_ENCODING creation option type from an input grib2 file for use in the output grib2 file

 

Thanks Shane/Sean. On the pygeoapi side, I also stumbled on this issue and was able to allow for user specific options in configuration and passing that on to rasterio accordingly. Definitely not dynamic, but worked for me, at least from the perspective of consistent (enough) production of GRIB2.

Cheers

..Tom


window_transform doesn't produce accurate coordinate transforms

pranay1117@...
 

Hi, I'm new to rasterio and working with geographical data. So far, it's been a great experience working with this library - but I've run into an issue when working with the Window module. My use case is as follows: 
- I obtain 2 bounding coordinates from the user (top-left and bottom-right) and convert them to (row, col) using `source.index(lon, lat)`.
- I create a window using these bounds and read the dataset from that. Also, I get the appropriate transform for this window.
selected_window = Window.from_slices((row0, row1), (col0, col1))
data = src.read(1, window=selected_window)
affine_transform = src.window_transform(selected_window)

Afterwards, I am using this `data` numpy array to select various points in these bounds and plot them on a map. However, the (lat, long) coordinates obtained using the affine_transform above are not very accurate. When I run the following:
print("Lat/Lon bounds: {}".format(src.window_bounds(selected_window)))
print("Transform at (0,0): {}".format(affine_transform * (0,0)))
r, c = data.shape
print("Transform at {}: {}".format((r-1, c-1), affine_transform * (r-1, c-1)))

The output is as follows:
Source lat/lon bounds: (72.45 20.69, 82.39, 28.78)
Transform at (0,0): (72.42, 28.78)
Transform at (970, 1193): (80.53, 18.841)

As is evident from the output, the latitude must lie between (20.69,28.78) but the transform value 18.841 at (970, 1193) is not in that interval. 

I saw some earlier posts in this group that suggested updating to version 1.1.5. I've done that but the issue persists. 
Any help will be greatly appreciated. Thanks!



rasterio opens file from AWS S3 bucket on local machine, but can't find file when deployed to Google App Engine

Judson Buescher
 

Hello,

 

I'm trying to run a python flask app, that opens a file hosted on an S3 bucket and analyzes the data. When I deploy the flask app on my local host, it runs like expected. Once I deploy the app to the cloud, it suddenly cannot find the file. I feel like this might be a bug, but i'm too afraid to post in the github. I'll attach the relevant snippet of code, and the error I'm receiving in the GCP Logs. I apologize for terrible formatting as I'm unaware of how to do it properly in this forum.
Here's the code: 

session = boto3.Session(aws_access_key_id=access_key, aws_secret_access_key=secret_access_key)    

def start_process(map, lon, lat):    
    with rasterio.Env(AWSSession(session)):        
        return getCoordinatePixel(map,lon,lat)  url = 's3://metcon-test-bucket/high_plus.tif'  

app = Flask(__name__) @app.route('/api')
def api_id():    
    lat = int(request.args['lat'])    
    lon = int(request.args['lon'])          
    return jsonify(start_process(url, lon, lat))  

if __name__ == '__main__':    
    app.run()

 and here's the error:



Re: Reprojecting GEOTIFF to epsg:3857 with scrambled CRS data (possibly GCJ-02)

Even Rouault
 

On mardi 21 juillet 2020 18:52:10 CEST toku2288 via groups.io wrote:

> Hi, Just an update. Turns GeoTIFF does use GCJ encoding, However, There seem

> to be some open-source projects online which can reverse the transformation

> to high accuracy, Is there some way to use this algorithm to reproject the

> GeoTIFF file?

 

This had been discussed for PROJ in

https://github.com/OSGeo/PROJ/issues/982

 

Summary: apparently technically doable, but for legal reasons for folks operating in China, cannot be included in PROJ.

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com


Reprojecting output image to input image crs where the image sizes differ

ashnair0007@...
 

I have a geotiff (size= 4490 x 7341) that I cut into 256 x 256 crops to do an operation. Once the operation is done, I merge them back together. Since the original dimensions aren't divisible by 256 the resulting image will be of size 4608 x 7424. However the output tif and the input tif are of the same location, so they need to be geo-referenced the same. Is it possible to do so?


Re: Reprojecting GEOTIFF to epsg:3857 with scrambled CRS data (possibly GCJ-02)

toku2288@...
 

Hi, Just an update. Turns GeoTIFF does use GCJ encoding, However, There seem to be some open-source projects online which can reverse the transformation to high accuracy, Is there some way to use this algorithm to reproject the GeoTIFF file?


Re: Convert BoundingBox data into lat/lon

toku2288@...
 

Hi Guillaume,

Thanks, It worked. I previously tried using the CRS EPSG:3857 using that method, and it returned lat/lon values out of range. Is there some way to get the EPSG:3857 values?


Re: Convert BoundingBox data into lat/lon

Guillaume Lostis
 

Hi,

The rasterio.warp.transform_bounds function should do what you want: https://rasterio.readthedocs.io/en/latest/api/rasterio.warp.html#rasterio.warp.transform_bounds

Pass it the CRS of your dataset as src_crs, CRS({'init': 'EPSG:4326'}) as dst_crs, and your bounds, and you will get your bounds in lon/lat as output.

Regards, 

Guillaume Lostis 

On Tue, 21 Jul 2020, 09:36 toku2288 via groups.io, <toku2288=protonmail.com@groups.io> wrote:
Hi, I need to convert the bounding box data of a TIF file to lat/lon co-ordinates so I can request data from another API

GridTIF.bounds returns a BoundingBox object like the following:
BoundingBox(left=-2638055.0614613006, bottom=340165.4847189281, right=2410944.9385386994, top=5922165.484718928)
Is there some way to convert this to lat/lon of the two corners in it's current CRS?
 


Convert BoundingBox data into lat/lon

toku2288@...
 

Hi, I need to convert the bounding box data of a TIF file to lat/lon co-ordinates so I can request data from another API

GridTIF.bounds returns a BoundingBox object like the following:
BoundingBox(left=-2638055.0614613006, bottom=340165.4847189281, right=2410944.9385386994, top=5922165.484718928)
Is there some way to convert this to lat/lon of the two corners in it's current CRS?
 


Reprojecting GEOTIFF to epsg:3857 with scrambled CRS data (possibly GCJ-02)

toku2288@...
 

Hello, I am trying to reproject a population distribution raster file of China into 3857 space. I tried but the images are not aligning. Looking at the CRS it has a lot of different standards, I am not sure what to make of it.

The meta of the source file looks something like this:
{'driver': 'GTiff', 
'dtype': 'int32', 'nodata': 2147483647.0, 'width': 5049, 'height': 5582, 'count': 1, 'crs': CRS.from_wkt('PROJCS["Albers_Conic_Equal_Area",GEOGCS["GCS_Unknown_datum_based_upon_the_Krassowsky_1940_ellipsoid",DATUM["Not_specified_based_on_Krassowsky_1940_ellipsoid",SPHEROID["Krasovsky_1940",6378245,298.3,AUTHORITY["EPSG","7024"]],AUTHORITY["EPSG","6024"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Albers_Conic_Equal_Area"],PARAMETER["standard_parallel_1",25],PARAMETER["standard_parallel_2",47],PARAMETER["latitude_of_center",0],PARAMETER["longitude_of_center",105],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]]]'), 'transform': Affine(1000.0, 0.0, -2638055.0614613006, 0.0, -1000.0, 5922165.484718928)}
PopTIF.crs.data returns:
{'proj': 'aea',
 'lat_1': 25,
 'lat_2': 47,
 'lat_0': 0,
 'lon_0': 105,
 'x_0': 0,
 'y_0': 0,
 'ellps': 'krass',
 'units': 'm',
 'no_defs': True}
after running the reprojection algorithm:
{'driver': 'GTiff',
 'dtype': 'int32',
 'nodata': 2147483647.0,
 'width': 6708,
 'height': 6102,
 'count': 1,
 'crs': CRS.from_dict(init='epsg:3857'),
 'transform': Affine(1184.7466003556895, 0.0, 7551104.353792087,
        0.0, -1184.7466003556895, 7299154.069454485)}




I am then using PopTIF.dataset_mask()  to mask the area of a Heightfield GEOTIFF (I requested the raster area of this file to be the same as the lat/lon of the PopTIF data:  'lat_1': 25,'lat_2': 47, 'lat_0': 0, 'lon_0': 105
 {'driver': 'GTiff', 
'dtype': 'float32', 'nodata': -32767.0, 'width': 2816, 'height': 2288, 'count': 1, 'crs': CRS.from_dict(init='epsg:4326'), 'transform': Affine(0.02197265625, 0.0, 73.4765625, 0.0, -0.02197265625, 53.7890625)}

 EPSG:4326 to EPSG:3857

{'driver': 'GTiff',
 'dtype': 'float32',
 'nodata': -32767.0,
 'width': 2594,
 'height': 2537,
 'count': 1,
 'crs': CRS.from_dict(init='epsg:3857'),
 'transform': Affine(2655.775407243914, 0.0, 8179373.5227401415,
        0.0, -2655.775407243914, 7130308.041557059)}


So, After moving both GEOTIFFs into 3857 space. I then needed to upscale the Heightfield so the array matched the the size of the Population data:

import skimage.transform as st

pathToRaster  = r'./DEMTIF/cnl10_3857.tif'
newsize = reprojectedmask.shape


with rio.open(pathToRaster) as src:
    arr = src.read()
    arr[arr ==-9999] = np.nan
    arr  = arr[0,:,:]
    new_shape = newsize
    newarr= st.resize(arr, new_shape, mode='constant')

Resulting array can be masked, and looks something like this, Notice how there is no height data in areas such as Taiwan, There has been a mismatch with the projections, However I cant tell where/how. This may be due to the fact I specified the lat/lon in AEA(plus other CRS)? Would appreciate any guidance on how to find the correct area to sample and how to transform it 


Re: Different values when I use a window

Sean Gillies
 

Hi,

We fixed a window shift bug in 1.1.5. See


I wonder if this is another bug or if somehow your installation didn't get the fix. Can you say what your source of the rasterio package is?

On Sun, Jul 5, 2020 at 5:11 PM <adrianocorbelinoii@...> wrote:
I am using rasterio 1.1.5.



--
Sean Gillies