Date   

rio-mbtiles 1.6a2

Sean Gillies
 

Hi all,

I've uploaded early pre-releases of rio-mbtiles 1.6 to PyPI. See https://github.com/mapbox/rio-mbtiles/blob/main/CHANGES.txt#L4 for detailed changes. In a nutshell, the new features are

* Support for 8-bit input and output
* Automatic alpha channel creation for RGB input when the --rgba option is used
* An option to keep empty tiles with no valid data

I'd love to hear reports about the pre-release.

--
Sean Gillies


Re: Creating WarpedVRT with multiple files

Sean Gillies
 

Hello,

On Thu, Jul 22, 2021 at 9:09 AM <remi.desgrange@...> wrote:
Hi,

It may be silly, but It does not seem possible to create a WarpedVRT object with multiple files ? Or is it ?

Thanks.

Rasterio's WarpedVRT class uses GDALCreateWarpedVRT (https://gdal.org/api/gdalwarp_cpp.html#_CPPv419GDALCreateWarpedVRT12GDALDatasetHiiPdP15GDALWarpOptions) and thus can wrap only a single dataset.

--
Sean Gillies


Creating WarpedVRT with multiple files

remi.desgrange@...
 

Hi,

It may be silly, but It does not seem possible to create a WarpedVRT object with multiple files ? Or is it ?

Thanks.


Re: No geoKeys generated when writing a geotiff

devsantiago
 

Hi,

We are currently using rasterio version 1.1.8 and it was installed following this procedure : https://rasterio.readthedocs.io/en/latest/installation.html#installing-with-anaconda


rio-color 1.0.4

Sean Gillies
 

Hi all,

In https://github.com/mapbox/rio-color/pull/79 it was pointed out that the wheels for 1.0.3 were built with versions of numpy that are too new for some users. In version 1.0.4 we depend on the "oldest-supported-numpy" package for Python versions >= 3.5. This means that, for example, the 1.0.4 wheels for Python 3.9 will be compatible with numpy versions >= 1.19. 

--
Sean Gillies


Understanding use of epsilon in transform.rowcol

tom.russell@...
 

Hi, I'm trying to understand the use of epsilon in rasterio.transform.rowcol

Specifically here - https://github.com/groutr/rasterio/blob/95d77c5e8553405a9dee373ae815072cf91d5cd1/rasterio/transform.py#L240 - where, as I read it, a point to be transformed from world coordinates to raster row/column indices is nudged by epsilon before being transformed, in the direction of the interior of the raster cell.

Why do we need to translate it that small amount before transforming? My guess is that it avoids later rounding or off-by-one errors, but I'm struggling to come up with an example or small test case.

Am I on the right lines? Is there good background material (or keywords to search on) that I'm missing for this kind of thing?

Thanks in advance, I hope this is the right forum for these questions.


Re: writing raster files produces black tifs

Julius Schmiedt
 

Thank you for the hint, I wasn´t aware of this ArcGIS behaviour concerning binary raster files. A first test was successful!

Am 12.07.2021 um 15:46 schrieb Hobart, Geordie (NRCan/RNCan):

Hi Julius.

 

If I understand you correctly, this is a Arc question.

IFF the values exist in your output file. You can test by loading your results file in a new process and rasterio.plot.show().

Arc tries to display a byte file for values between 0-255, hence no contrast between 0 and 1. In Arc,  open the layer properties of your  new GTiff and adjust the display values under the symbology tab.  

 

If not and the results file is all zeros then it is probably and error in your code. I had a similar issue where I incorrectly set the value NBITS =1 when calling rasterio.open() and on linux systems my output was empty.

 

Good luck.

Geordie

 

 

From: main@rasterio.groups.io <main@rasterio.groups.io> On Behalf Of Julius Schmiedt
Sent: July 12, 2021 02:18
To: main@rasterio.groups.io
Subject: [rasterio] writing raster files produces black tifs

 

Hello everybody,

I´m quite new to rasterio and Python, but my task is, to produce binary raster files based on the information of sentinel2 satellite images (GTiff). Therefore I use Spyder IDE and the numpy.where() method to calculate new arrays and save them as GTiffs.

Although the output arrays contain 0 and 1 values, the output tif only shows black (in the shape of the original gtiff), if opened with standard image viewer or ArcGIS Desktop. Only if opened via Python script in Spyder, I am able to display the output tifs as expected (using rasterio.plot.show() method).

If anybody has an idea, what the problem might be, I would appreciate any help.

Thanks!

--
Julius Schmiedt
Wissenschaftliche Hilfskraft

Department Landschaftsökologie
Helmholtz-Zentrum für Umweltforschung GmbH - UFZ
Permoserstraße 15 | 04318 Leipzig
Tel.: +49 341 235 - 1058
Email: julius.schmiedt@...
Web: https://www.ufz.de

--
Julius Schmiedt
Wissenschaftliche Hilfskraft

Department Landschaftsökologie
Helmholtz-Zentrum für Umweltforschung GmbH - UFZ
Permoserstraße 15 | 04318 Leipzig
Tel.: +49 341 235 - 1058
Email: julius.schmiedt@...
Web: https://www.ufz.de


Re: No geoKeys generated when writing a geotiff

Sean Gillies
 

Hi,

On Mon, Jul 12, 2021 at 7:36 AM devsantiago <santiago@...> wrote:
Could be, but i doubt it. The other geotiffs I'm able to successfully read (geoKeys) and display inside the web application (using npm's geotiff package) are generated using professional GIS software. Further question, in my full script, there are other places (not shown in sample code) where a dataset is open with rasterio and some information is read to calculate the overlay. Could it be that the geotiffs that i'm using to generate my output image don't have any information encoded, which would explain why nothing is passed to my output image? Also, is there some sample code somewhere using rasterio to generate a geotiff (and its geoKeys) that I could look at ?
,_._,_

I'm stumped. Rasterio doesn't have any native understanding of GeoTIFF geokeys. To define a dataset's georeferencing, we call GDALSetGeoTransform and GDALSetProjection and GDAL's GeoTIFF driver and libgeotiff take care of the details. We haven't had a report of rasterio (and GDAL) writing a non-conforming GeoTIFF in a long time. The tutorial at https://rasterio.readthedocs.io/en/latest/quickstart.html#opening-a-dataset-in-writing-mode produces a GeoTIFF that can be used with GDAL and QGIS. Can you tell me what rasterio version you are using and what its provenance is? PyPI? Anaconda? Something else?

--
Sean Gillies


Re: writing raster files produces black tifs

Hobart, Geordie (NRCan/RNCan)
 

Hi Julius.

 

If I understand you correctly, this is a Arc question.

IFF the values exist in your output file. You can test by loading your results file in a new process and rasterio.plot.show().

Arc tries to display a byte file for values between 0-255, hence no contrast between 0 and 1. In Arc,  open the layer properties of your  new GTiff and adjust the display values under the symbology tab.  

 

If not and the results file is all zeros then it is probably and error in your code. I had a similar issue where I incorrectly set the value NBITS =1 when calling rasterio.open() and on linux systems my output was empty.

 

Good luck.

Geordie

 

 

From: main@rasterio.groups.io <main@rasterio.groups.io> On Behalf Of Julius Schmiedt
Sent: July 12, 2021 02:18
To: main@rasterio.groups.io
Subject: [rasterio] writing raster files produces black tifs

 

Hello everybody,

I´m quite new to rasterio and Python, but my task is, to produce binary raster files based on the information of sentinel2 satellite images (GTiff). Therefore I use Spyder IDE and the numpy.where() method to calculate new arrays and save them as GTiffs.

Although the output arrays contain 0 and 1 values, the output tif only shows black (in the shape of the original gtiff), if opened with standard image viewer or ArcGIS Desktop. Only if opened via Python script in Spyder, I am able to display the output tifs as expected (using rasterio.plot.show() method).

If anybody has an idea, what the problem might be, I would appreciate any help.

Thanks!

--
Julius Schmiedt
Wissenschaftliche Hilfskraft

Department Landschaftsökologie
Helmholtz-Zentrum für Umweltforschung GmbH - UFZ
Permoserstraße 15 | 04318 Leipzig
Tel.: +49 341 235 - 1058
Email: julius.schmiedt@...
Web: https://www.ufz.de


Re: No geoKeys generated when writing a geotiff

devsantiago
 

Could be, but i doubt it. The other geotiffs I'm able to successfully read (geoKeys) and display inside the web application (using npm's geotiff package) are generated using professional GIS software. Further question, in my full script, there are other places (not shown in sample code) where a dataset is open with rasterio and some information is read to calculate the overlay. Could it be that the geotiffs that i'm using to generate my output image don't have any information encoded, which would explain why nothing is passed to my output image? Also, is there some sample code somewhere using rasterio to generate a geotiff (and its geoKeys) that I could look at ?


writing raster files produces black tifs

Julius Schmiedt
 

Hello everybody,

I´m quite new to rasterio and Python, but my task is, to produce binary raster files based on the information of sentinel2 satellite images (GTiff). Therefore I use Spyder IDE and the numpy.where() method to calculate new arrays and save them as GTiffs.

Although the output arrays contain 0 and 1 values, the output tif only shows black (in the shape of the original gtiff), if opened with standard image viewer or ArcGIS Desktop. Only if opened via Python script in Spyder, I am able to display the output tifs as expected (using rasterio.plot.show() method).

If anybody has an idea, what the problem might be, I would appreciate any help.

Thanks!

--
Julius Schmiedt
Wissenschaftliche Hilfskraft

Department Landschaftsökologie
Helmholtz-Zentrum für Umweltforschung GmbH - UFZ
Permoserstraße 15 | 04318 Leipzig
Tel.: +49 341 235 - 1058
Email: julius.schmiedt@...
Web: https://www.ufz.de


Re: MODIS (HDF4) data with MemoryFile

Hussain Hassam
 

Ahh I see. Thank you Sean & Evan for your quick responses. We will proceed with our backup plan of storing data to disk or bucket and using it from there!

Best Regards,
Hussain

-----Original Message-----
From: main@rasterio.groups.io <main@rasterio.groups.io> On Behalf Of Even Rouault via groups.io
Sent: July 9, 2021 3:49 PM
To: main@rasterio.groups.io
Subject: Re: [rasterio] MODIS (HDF4) data with MemoryFile


It may be possible that in-memory HDF4 is not supported by GDAL.
in-memory HDF4 or any /vsi virtual file system will not work with HDF4.
The underlying library has no abstraction for I/O and uses plain C standard library API


--
https://can01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.spatialys.com%2F&;data=04%7C01%7Chhassam%40hatfieldgroup.com%7C266bcd2bcf6547425a5808d9432bb181%7Ce0b2a496c1864e92b4b07d466b05d8d9%7C0%7C0%7C637614677490906436%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000&amp;sdata=Mzer7qw088AOX0lDrFS3R3P%2BWk8an9y9cs8WGDHrDvM%3D&amp;reserved=0
My software is free, but my time generally not.


Re: MODIS (HDF4) data with MemoryFile

Even Rouault
 

It may be possible that in-memory HDF4 is not supported by GDAL.
in-memory HDF4 or any /vsi virtual file system will not work with HDF4. The underlying library has no abstraction for I/O and uses plain C standard library API


--
http://www.spatialys.com
My software is free, but my time generally not.


Re: MODIS (HDF4) data with MemoryFile

Sean Gillies
 

Hi Hussain,

Thanks for posting here. I don't participate in that Slack site and so I don't see any rasterio questions there.

Can you show the exact code you are using?

I see in GDAL's tests and in an announcement https://lists.osgeo.org/pipermail/gdal-dev/2018-August/048934.html that in-memory HDF5 is supported. It may be possible that in-memory HDF4 is not supported by GDAL. I am not a regular HDF user and so I am not certain whether that should be expected to work or not.

On Fri, Jul 9, 2021 at 4:29 PM Hussain Hassam via groups.io <hhassam=hatfieldgroup.com@groups.io> wrote:

Hello,

 

I posted the below message in the Spatial Community #mapbox slack channel, but I have yet to receive a response. I am hoping to get some information about this before I consider it a bug and post on github issues.

 

We have been trying to open HDF4 files (MODIS data) directly from a source (e.g. https://e4ftl01.cr.usgs.gov//DP131/MOLT/MOD09A1.061/2021.06.26/MOD09A1.A2021177.h14v00.061.2021186035450.hdf) with rasterio (in Python).The first issue is that the source requires authentication. We solved by using requests to download the HDF file into a buffer (python io.BytesIO object. We also seeked it back to 0,0 after writing to it). Our thought was to then load this into rasterio using rasterio.io.MemoryFile. When we try to open the memoryfile, we get an error:

 

CPLE_OpenFailedError: '/vsimem/c604e9d5-cc4b-42bc-a34a-de4e541355cc/c604e9d5-cc4b-42bc-a34a-de4e541355cc.hdf' not recognized as a supported file format.

 

The odd part about this is if we write the io.BytesIO object into a file on disk first (example.hdf) then do a rasterio.open('example.hdf'), everything works as expected, however, we are trying to avoid any sort of IO for a streamlined process.

 

Thank you for your time,

Hussain Hassam

_._,_._,_

--
Sean Gillies


MODIS (HDF4) data with MemoryFile

Hussain Hassam
 

Hello,

 

I posted the below message in the Spatial Community #mapbox slack channel, but I have yet to receive a response. I am hoping to get some information about this before I consider it a bug and post on github issues.

 

We have been trying to open HDF4 files (MODIS data) directly from a source (e.g. https://e4ftl01.cr.usgs.gov//DP131/MOLT/MOD09A1.061/2021.06.26/MOD09A1.A2021177.h14v00.061.2021186035450.hdf) with rasterio (in Python).The first issue is that the source requires authentication. We solved by using requests to download the HDF file into a buffer (python io.BytesIO object. We also seeked it back to 0,0 after writing to it). Our thought was to then load this into rasterio using rasterio.io.MemoryFile. When we try to open the memoryfile, we get an error:

 

CPLE_OpenFailedError: '/vsimem/c604e9d5-cc4b-42bc-a34a-de4e541355cc/c604e9d5-cc4b-42bc-a34a-de4e541355cc.hdf' not recognized as a supported file format.

 

The odd part about this is if we write the io.BytesIO object into a file on disk first (example.hdf) then do a rasterio.open('example.hdf'), everything works as expected, however, we are trying to avoid any sort of IO for a streamlined process.

 

Thank you for your time,

Hussain Hassam


Re: No geoKeys generated when writing a geotiff

Sean Gillies
 

Hi,

Right, the typo wouldn't make a difference since "STANDARD" is already the GDAL GeoTIFF default.

At https://github.com/geotiffjs/geotiff.js I see "Rudimentary extraction of geospatial metadata". Perhaps the package isn't complete in this area? I didn't immediately see tests of reading geotiff keys.

On Thu, Jul 8, 2021 at 10:34 AM devsantiago <santiago@...> wrote:

Hi Sean, 

Thanks for your quick answer and for noting this, i had indeed a typo. However, this did not solve my problem. Am I missing something else ? There is not much info in the documentation and I feel like I included everything needed.

._,_._,

--
Sean Gillies


Re: No geoKeys generated when writing a geotiff

devsantiago
 

Hi Sean, 

Thanks for your quick answer and for noting this, i had indeed a typo. However, this did not solve my problem. Am I missing something else ? There is not much info in the documentation and I feel like I included everything needed.


Re: Open EHdr raster with hdr header file from memory

fkluibenschaedl@...
 

On Fri, Jul 2, 2021 at 01:12 AM, Luke wrote:
Replacing your FTP code with direct access using the /vsicurl/ virtual filesystem:

This is a very elegant solution, thanks for pointing it out and thanks for the sample code! The only drawback I see is that it is not possible to access multiple products (files within the archive) while at the same time downloading the archive only once.

Thanks again,
F


WarpedVRT and mixed data types

Thomas Maschler
 

Hi,

I noticed some odd behavior when using WarpedVRT in combination with a VRT with mixed Data Types. 

When trying to open the mixed data type raster using rasterio, it throws an error due to no support for mixed data types. However, when trying to open the same dataset as WarpedVRT, I do not get such error. Instead it returns an all-zero array.

I would either expect to get the same error message, receive an array with all data cast to the largest data type or receive an array with mixed data types.

 

import rasterio
import subprocess
from affine import Affine

from rasterio.vrt import WarpedVRT



# First raster saved as uint16
a = np.array([[0,1,2,3,4]]).astype("uint16")

with rasterio.open("uint16.tif", "w", width=5, height=1, count=1, dtype=rasterio.uint16, nodata=0) as dst:

    dst.write(a, 1)

# Second raster saved as int32
b = np.array([[1,2,3,4,0]]).astype("int32")

with rasterio.open("int32.tif", "w", width=5, height=1, count=1, dtype=rasterio.int32, nodata=0) as dst:

    dst.write(b, 1)

# create a VRT combining both datasets into a two band raster using gdalbuildvrt

subprocess.run(["gdalbuildvrt", "-separate", "mixed_datatype.tif", "uint16.tif", "int32.tif"], capture_output=True)

# Try to open the two band raster using rasterio

with rasterio.open("mixed_datatype.tif") as src:

    print(src.profile)

    print(src.read())  # This throws an error due to mixed data types

{'driver': 'VRT', 'dtype': 'uint16', 'nodata': 0.0, 'width': 5, 'height': 1, 'count': 2, 'crs': None, 'transform': Affine(1.0, 0.0, 0.0, 0.0, 1.0, 0.0), 'tiled': False}

Traceback (most recent call last):

  File "<input>", line 3, in <module>

  File "rasterio/_io.pyx", line 278, in rasterio._io.DatasetReaderBase.read

ValueError: more than one 'dtype' found

 

# Try to open the two band raster using WarpedVRT

with rasterio.open("mixed_datatype.tif") as src:

    with WarpedVRT(src, crs="EPSG:4326", transform=transform, width=5, height=1) as vrt:

        print(vrt.profile)

        print(vrt.read())  # This does not throw an error but returns an array with only zeros

    

{'driver': 'VRT', 'dtype': 'int32', 'nodata': 0.0, 'width': 5, 'height': 1, 'count': 2, 'crs': CRS.from_epsg(4326), 'transform': Affine(1.0, 0.0, 0.0,

       0.0, -1.0, 0.0), 'tiled': False}

[[[0 0 0 0 0]]

 [[0 0 0 0 0]]]




Re: How are calculate_default_transform and aligned_target different?

amaury.dehecq@...
 

Hi Sean,

Ok thanks for the clarifications. I think I found GDAL's logic around here: https://github.com/OSGeo/gdal/blob/d01976cd4fc11144662f44ec4e658ce1871f1713/gdal/apps/gdalwarp_lib.cpp#L3105

It would indeed be useful to have a function that returns the best transform and shape given any bounds, resolution and shape. That's what I tried to do here in a Python implementation: https://github.com/GlacioHack/GeoUtils/blob/12bd607fe8024720ff264545abf0cc5ba0156ea6/geoutils/georaster.py#L879 (until L 932), if that can be of any help. In the meantime, I think I'll just keep my current implementation, which seems to provide similar results to what gdalwarp would do (I'll check the code in detail later).

Amaury

 

1 - 20 of 828