Date   
Mask to alpha band

Denis Rykov
 

Hi, folks!

Is it possible to convert internal mask to alpha band? 

Re: Reproject with mask preserving

Sean Gillies
 

Hi Denis,

As you discovered in your post to gdal-dev, GDAL's warper doesn't output internal (or sidecar) mask bands. You must materialize these into an alpha band and then copy the data to produce a new internal mask band.

Here's a test that demonstrates how to create an internal mask band:


While putting together an example of calling reproject() with alpha band creation, I found a bug


When Rasterio 1.0.1 is published to PyPI, you will be able to call reproject() and materialize source's mask band into the destination alpha band specified by the dst_alpha keyword argument as shown in this test:



On Sat, Jul 21, 2018 at 4:45 PM, Denis Rykov <rykovd@...> wrote:
Hi there!

I have GeoTIFF files with internal mask and I need to reproject them and not lose mask.

I've found that direct generating warped mask band is not supported by GDAL: https://github.com/OSGeo/gdal/issues/689

rasterio.warp.reproject() also doesn't preserve mask.

How I can use rasterio for achieving what I need? I'm looking for something similar described in this issue but with rasterio.




--
Sean Gillies

Rasterio 1.0.1

Sean Gillies
 

Hi all,

Rasterio 1.0.1 is on PyPI now. It fixes 4 bugs found since 1.0.0:

- Bounding envelopes are densified at higher precision in transform_bounds to
  fix #1411.
- Add LERC compression to enums.Compression (#1412).
- The reproject() function now passes dst_alpha properly to _reproject(), which
  unlocks materialization of warp destination alpha band masks (#1417).
- The --dimensions and --src-bounds options of rio-warp can be used together
  as expected (#1418).

--
Sean Gillies

Re: Mask to alpha band

Sean Gillies
 

Hi Denis,

There's no way to do it in Rasterio without creating a new dataset and copying the band and mask arrays from the original.

On Mon, Jul 23, 2018 at 2:08 PM, Denis Rykov <rykovd@...> wrote:
Hi, folks!

Is it possible to convert internal mask to alpha band? 




--
Sean Gillies

read masked and overviews when no-data value is defined

Guy Doulberg
 

Hi all,

We are trying to read masked array using `DatasetReaderBase.read`, reading the documentation it looks that if we pass masked=True we will get what we wanted. 

I was testing the read on raster with pre-built overviews and I used the overviews in read. From functionality point of view, I always got the right result, but I was tesing the data block that are being used, and I noticed this behaviors:


1. When using band-mask, there was a read only from the overview to construct the mask -- as excpected
2. When using no-data value there were reads the data of the original resolution and not from the nearest overview? when mask=False it reads data only from the overview,  why is that?

A workaround for that is to read with no mask and to construct the mask myself using the no-data value.

Thanks.

`import rasterio` error on Win7 and Mac: ValueError: numpy.dtype has the wrong size, try recompiling

nurbachris@...
 

I am having the same problem on Windows 7 and Mac OS X High Sierra:

1) Windows 7 (64 bits)

Using python 2.7 (32 bits).

I have installed GDAL Core and its Python bindings using the executables here (resp. gdal-202-1500-core.msi and GDAL-2.2.3.win32-py2.7.msi). This worked well and I could successfully use gdal in python scripts.

Then I tried to install rasterio using rasterio‑1.0.3‑cp27‑cp27m‑win32.whl:
py -2.7 -m pip install rasterio‑1.0.3‑cp27‑cp27m‑win32.whl

This went through successfully.

However when starting python and invoking import rasterio I get the following error message:

Traceback (most recent call last):
  File "test1.py", line 1, in <module>
    import rasterio
  File "C:\Python27\lib\site-packages\rasterio\__init__.py", line 31, in <module>
    from rasterio._base import gdal_version
  File "rasterio\_base.pyx", line 1, in init rasterio._base
  File "__init__.pxd", line 164, in init rasterio._shim
ValueError: numpy.dtype has the wrong size, try recompiling. Expected 52, got 56

2) Mac OS X High Sierra (10.13.6)

Interestingly, I am having exactly the same error on Mac OS X High Sierra after installing GDAL Complete 2.1 from here, and rasterio using:
pip2.7 install rasterio.

The installation steps are successful, but I get the same error as above when executing import rasterio in python 2.7:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/Library/Python/2.7/site-packages/rasterio/__init__.py", line 23, in <module>
    from rasterio._base import gdal_version
  File "rasterio/_base.pyx", line 1, in init rasterio._base
  File "__init__.pxd", line 164, in init rasterio._shim
ValueError: numpy.dtype has the wrong size, try recompiling. Expected 88, got 96

I have read posts advising to upgrade numpy, however it is installed as distutils project (on both systems) and therefore I cannot upgrade with pip. So I am not sure what is the best way forward.

Any help would be greatly appreciated.

Re: `import rasterio` error on Win7 and Mac: ValueError: numpy.dtype has the wrong size, try recompiling

Sean Gillies
 

I am sorry, but upgrading numpy is the only way forward. If you don't want to overwrite your existing numpy installation, you have several  options. I believe the following two are the best. Other group members may have different experiences.

1. Download python 3.6 or 3.7 from python.org and install it. Now you have a brand new (and better) version of python into which you can install packages. On your Mac, `python3 -m pip install rasterio` will cause numpy to be fetched and installed, too. On your Windows computer, you will need to get numpy separately from https://www.lfd.uci.edu/~gohlke/pythonlibs/#numpy and install it.

2. Install Miniconda (3.6 again) https://conda.io/miniconda.html and get Rasterio from conda-forge: `conda install -c conda-forge rasterio`.

Good luck!

On Tue, Aug 7, 2018 at 8:33 AM <nurbachris@...> wrote:

I am having the same problem on Windows 7 and Mac OS X High Sierra:

1) Windows 7 (64 bits)

Using python 2.7 (32 bits).

I have installed GDAL Core and its Python bindings using the executables here (resp. gdal-202-1500-core.msi and GDAL-2.2.3.win32-py2.7.msi). This worked well and I could successfully use gdal in python scripts.

Then I tried to install rasterio using rasterio‑1.0.3‑cp27‑cp27m‑win32.whl:
py -2.7 -m pip install rasterio‑1.0.3‑cp27‑cp27m‑win32.whl

This went through successfully.

However when starting python and invoking import rasterio I get the following error message:

Traceback (most recent call last):
  File "test1.py", line 1, in <module>
    import rasterio
  File "C:\Python27\lib\site-packages\rasterio\__init__.py", line 31, in <module>
    from rasterio._base import gdal_version
  File "rasterio\_base.pyx", line 1, in init rasterio._base
  File "__init__.pxd", line 164, in init rasterio._shim
ValueError: numpy.dtype has the wrong size, try recompiling. Expected 52, got 56

2) Mac OS X High Sierra (10.13.6)

Interestingly, I am having exactly the same error on Mac OS X High Sierra after installing GDAL Complete 2.1 from here, and rasterio using:
pip2.7 install rasterio.

The installation steps are successful, but I get the same error as above when executing import rasterio in python 2.7:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/Library/Python/2.7/site-packages/rasterio/__init__.py", line 23, in <module>
    from rasterio._base import gdal_version
  File "rasterio/_base.pyx", line 1, in init rasterio._base
  File "__init__.pxd", line 164, in init rasterio._shim
ValueError: numpy.dtype has the wrong size, try recompiling. Expected 88, got 96

I have read posts advising to upgrade numpy, however it is installed as distutils project (on both systems) and therefore I cannot upgrade with pip. So I am not sure what is the best way forward.

Any help would be greatly appreciated.



--
Sean Gillies

Re: read masked and overviews when no-data value is defined

Sean Gillies
 

Hi Guy,

A month ago we resolved an issue in Rasterio that has some commonalities with yours. There was a report that overviews were being missed when reading masked arrays and it turned out that the overviews themselves were in error.


In Rasterio, we will skip reading a mask from the dataset if the dataset has no mask, but otherwise we don't call GDAL functions differently in the band mask vs nodata cases.

Can you tell us what nodata value is implicated in this problem and what software produced the data files?


On Mon, Aug 6, 2018 at 8:52 AM <guyd@...> wrote:
Hi all,

We are trying to read masked array using `DatasetReaderBase.read`, reading the documentation it looks that if we pass masked=True we will get what we wanted. 

I was testing the read on raster with pre-built overviews and I used the overviews in read. From functionality point of view, I always got the right result, but I was tesing the data block that are being used, and I noticed this behaviors:


1. When using band-mask, there was a read only from the overview to construct the mask -- as excpected
2. When using no-data value there were reads the data of the original resolution and not from the nearest overview? when mask=False it reads data only from the overview,  why is that?

A workaround for that is to read with no mask and to construct the mask myself using the no-data value.

Thanks.



--
Sean Gillies

Re: read masked and overviews when no-data value is defined

Guy Doulberg
 

He Sean,

We are creating  cog files using our code that takes as input any raster, tiles it add overviews and rearrange the tif to be a cog.

For doing so we are using rasetio shutils copy and rasterio build_overview.

The source rasters that we are working on have nodata

We did a variant of code in which we prior for creating the overviews we also created Band Mask

1. when we had no-data values, and read using the overviews it looked like it was constructing the overview from the original size

On Wed, Aug 8, 2018 at 12:14 AM, Sean Gillies <sean.gillies@...> wrote:
Hi Guy,

A month ago we resolved an issue in Rasterio that has some commonalities with yours. There was a report that overviews were being missed when reading masked arrays and it turned out that the overviews themselves were in error.


In Rasterio, we will skip reading a mask from the dataset if the dataset has no mask, but otherwise we don't call GDAL functions differently in the band mask vs nodata cases.

Can you tell us what nodata value is implicated in this problem and what software produced the data files?


On Mon, Aug 6, 2018 at 8:52 AM <guyd@...> wrote:
Hi all,

We are trying to read masked array using `DatasetReaderBase.read`, reading the documentation it looks that if we pass masked=True we will get what we wanted. 

I was testing the read on raster with pre-built overviews and I used the overviews in read. From functionality point of view, I always got the right result, but I was tesing the data block that are being used, and I noticed this behaviors:


1. When using band-mask, there was a read only from the overview to construct the mask -- as excpected
2. When using no-data value there were reads the data of the original resolution and not from the nearest overview? when mask=False it reads data only from the overview,  why is that?

A workaround for that is to read with no mask and to construct the mask myself using the no-data value.

Thanks.


--
Sean Gillies


Re: read masked and overviews when no-data value is defined

Guy Doulberg
 

clicked send too soon:

2. when using masked overviews, it looked like it was using the overview


later today I will try to run tests again

Do you care on how we created the original raster?


On Wed, Aug 8, 2018 at 11:43 AM, Guy Doulberg <guyd@...> wrote:
He Sean,

We are creating  cog files using our code that takes as input any raster, tiles it add overviews and rearrange the tif to be a cog.

For doing so we are using rasetio shutils copy and rasterio build_overview.

The source rasters that we are working on have nodata

We did a variant of code in which we prior for creating the overviews we also created Band Mask

1. when we had no-data values, and read using the overviews it looked like it was constructing the overview from the original size

On Wed, Aug 8, 2018 at 12:14 AM, Sean Gillies <sean.gillies@...> wrote:
Hi Guy,

A month ago we resolved an issue in Rasterio that has some commonalities with yours. There was a report that overviews were being missed when reading masked arrays and it turned out that the overviews themselves were in error.


In Rasterio, we will skip reading a mask from the dataset if the dataset has no mask, but otherwise we don't call GDAL functions differently in the band mask vs nodata cases.

Can you tell us what nodata value is implicated in this problem and what software produced the data files?


On Mon, Aug 6, 2018 at 8:52 AM <guyd@...> wrote:
Hi all,

We are trying to read masked array using `DatasetReaderBase.read`, reading the documentation it looks that if we pass masked=True we will get what we wanted. 

I was testing the read on raster with pre-built overviews and I used the overviews in read. From functionality point of view, I always got the right result, but I was tesing the data block that are being used, and I noticed this behaviors:


1. When using band-mask, there was a read only from the overview to construct the mask -- as excpected
2. When using no-data value there were reads the data of the original resolution and not from the nearest overview? when mask=False it reads data only from the overview,  why is that?

A workaround for that is to read with no mask and to construct the mask myself using the no-data value.

Thanks.


--
Sean Gillies



rio stack with compression

Luke Pinner
 

I stacked some singleband rasters into a multiband raster using `rio stack` with a compress=deflate creation option.  The compressed output filesize was approx 1.4x the sum of the uncompressed input filesizes.  Specifying --co interleave=band when running `rio stack` got the output filesize down to 0.25x the uncompressed input (per issue #70).   Running `gdalbuildvrt` to stack them, then `rio convert` or `gdal_translate` on the VRT with the same creation options as the original rio stack command also results in an output raster of around 0.25x the uncompressed input, but with default pixel interleaving.

I realise you're planning to document issues with compression (#77) but in the meantime, do you have any idea why compression is so poor with default pixel interleaving using `rio stack` but quite ok when using `rio convert` (with default pixel interleaving) on an already stacked VRT?

Luke

Re: rio stack with compression

Sean Gillies
 

Hi Luke,

I don't understand the issue entirely, but as Even Rouault explains in https://lists.osgeo.org/pipermail/gdal-dev/2015-June/041917.html a "naive" use of the write() method of a Rasterio dataset (and rio-stack does use it naively) can result in redundant blocks of data written to the GeoTIFF file. The workarounds are to specify band interleaving (as you discovered), or increase GDAL_CACHEMAX to allow the entire output file to be kept in memory until it is written. Remember: a dataset's write() method causes data to be written to the GDAL block cache and data is written from that cache to the GeoTIFF file when the block is bumped out by new data coming in, or when the dataset object is closed.

These are the issues that I need to bring up in Rasterio's documentation.

Hope this helps!

On Thu, Sep 6, 2018 at 5:37 AM Luke Pinner <lukepinnerau@...> wrote:
I stacked some singleband rasters into a multiband raster using `rio stack` with a compress=deflate creation option.  The compressed output filesize was approx 1.4x the sum of the uncompressed input filesizes.  Specifying --co interleave=band when running `rio stack` got the output filesize down to 0.25x the uncompressed input (per issue #70).   Running `gdalbuildvrt` to stack them, then `rio convert` or `gdal_translate` on the VRT with the same creation options as the original rio stack command also results in an output raster of around 0.25x the uncompressed input, but with default pixel interleaving.

I realise you're planning to document issues with compression (#77) but in the meantime, do you have any idea why compression is so poor with default pixel interleaving using `rio stack` but quite ok when using `rio convert` (with default pixel interleaving) on an already stacked VRT?

Luke



--
Sean Gillies

Rasterio 1.0.3

Sean Gillies
 

Hi all,

Rasterio 1.0.3 is on PyPI now, finally. It fixes 2 bugs found since 1.0.2:

- A regression in GeoJSON source handling in rio-rasterize (#1425) has been fixed.
- Rasterization of linear rings (#1431) is now supported.

Furthermore, we've made a switch to a multibuild-based infrastructure for building wheels: https://github.com/sgillies/rasterio-wheels. Please let me know if the new wheels don't behave as expected.

--
Sean Gillies

Rasterio 1.0.3.post1

Sean Gillies
 

Hi all,

The 1.0.3 wheels I built for OSX and uploaded to PyPI included a GDAL library with no curl support. I've made a post-release that fixes this problem. If you installed one of these wheels, please `pip install -U pip` to get 1.0.3.post1.

--
Sean Gillies

Re: Rasterio 1.0.3.post1

Denis Rykov
 

Hi Sean!

After updating rasterio I get the following error:

In [1]: from rasterio import open
 
In [2]: open('tests/data/raster/rgb.jp2')
---------------------------------------------------------------------------
CPLE_OpenFailedError                      Traceback (most recent call last)
rasterio/_base.pyx in rasterio._base.DatasetBase.__init__()
 
rasterio/_shim.pyx in rasterio._shim.open_dataset()
 
rasterio/_err.pyx in rasterio._err.exc_wrap_pointer()
 
CPLE_OpenFailedError: 'tests/data/raster/rgb.jp2' not recognized as a supported file format.
 
During handling of the above exception, another exception occurred:
 
RasterioIOError                           Traceback (most recent call last)
<ipython-input-2-19e7b16f33e4> in <module>()
----> 1 open('tests/data/raster/rgb.jp2')
 
~/sandbox/telluric/env/lib/python3.5/site-packages/rasterio/env.py in wrapper(*args, **kwds)
    401             else:
    402                 pass
--> 403             return f(*args, **kwds)
    404 
    405     return wrapper
 
~/sandbox/telluric/env/lib/python3.5/site-packages/rasterio/__init__.py in open(fp, mode, driver, width, height, count, crs, transform, dtype, nodata, sharing, **kwargs)
    215         # None.
    216         if mode == 'r':
--> 217             s = DatasetReader(path, driver=driver, **kwargs)
    218         elif mode == 'r+':
    219             s = get_writer_for_path(path)(path, mode, driver=driver, **kwargs)
 
rasterio/_base.pyx in rasterio._base.DatasetBase.__init__()
 
RasterioIOError: 'tests/data/raster/rgb.jp2' not recognized as a supported file format.

All works without any problems with rasterio 1.0.2:

In [1]: from rasterio import open
 
In [2]: open('tests/data/raster/rgb.jp2')
Out[2]: <open DatasetReader name='tests/data/raster/rgb.jp2' mode='r'>

Re: Rasterio 1.0.3.post1

Sean Gillies
 

Thanks for the report, Denis. I had upgraded GDAL to 2.3.1 in these wheel builds, and was blaming my inability to access these jp2 files on a lack of credentials, but maybe something has changed between Rasterio and JP2 in 2.3.0? I'm going to revert to 2.2.4, the version in the 1.0.2 wheels.


On Mon, Sep 10, 2018 at 7:10 AM Denis Rykov <rykovd@...> wrote:
Hi Sean!

After updating rasterio I get the following error:

In [1]: from rasterio import open
 
In [2]: open('tests/data/raster/rgb.jp2')
---------------------------------------------------------------------------
CPLE_OpenFailedError                      Traceback (most recent call last)
rasterio/_base.pyx in rasterio._base.DatasetBase.__init__()
 
rasterio/_shim.pyx in rasterio._shim.open_dataset()
 
rasterio/_err.pyx in rasterio._err.exc_wrap_pointer()
 
CPLE_OpenFailedError: 'tests/data/raster/rgb.jp2' not recognized as a supported file format.
 
During handling of the above exception, another exception occurred:
 
RasterioIOError                           Traceback (most recent call last)
<ipython-input-2-19e7b16f33e4> in <module>()
----> 1 open('tests/data/raster/rgb.jp2')
 
~/sandbox/telluric/env/lib/python3.5/site-packages/rasterio/env.py in wrapper(*args, **kwds)
    401             else:
    402                 pass
--> 403             return f(*args, **kwds)
    404 
    405     return wrapper
 
~/sandbox/telluric/env/lib/python3.5/site-packages/rasterio/__init__.py in open(fp, mode, driver, width, height, count, crs, transform, dtype, nodata, sharing, **kwargs)
    215         # None.
    216         if mode == 'r':
--> 217             s = DatasetReader(path, driver=driver, **kwargs)
    218         elif mode == 'r+':
    219             s = get_writer_for_path(path)(path, mode, driver=driver, **kwargs)
 
rasterio/_base.pyx in rasterio._base.DatasetBase.__init__()
 
RasterioIOError: 'tests/data/raster/rgb.jp2' not recognized as a supported file format.

All works without any problems with rasterio 1.0.2:

In [1]: from rasterio import open
 
In [2]: open('tests/data/raster/rgb.jp2')
Out[2]: <open DatasetReader name='tests/data/raster/rgb.jp2' mode='r'>



--
Sean Gillies

Fast X-array rasterization

Philipp Kats
 

Hi! I am trying to produce a map using rasterio's rasterization and datashader, following this gist

It works perfectly well, when I rasterize all the features at once, into one 2d numpy array;
However, what I want now is to produce a "cube", having separate raster for each feature separately, with a corresponding value (say, number of category). Plan is then to perform aggregation - say, take moda, from this cube.

However, that takes a lot of time on my machine (~3 hours vs 30 sec for one raster map), so I wonder if there are any suggestions/hacks on how I can speed it up.

Thank you in advance! 



Re: Fast X-array rasterization

Brendan Ward
 

Hi Philipp,

How are you looping over each value to create a new raster?  Are you able to share a gist of what you are doing?

Re: Rasterio 1.0.3.post1

Even Rouault
 

On lundi 10 septembre 2018 14:29:39 CEST Sean Gillies wrote:
Thanks for the report, Denis. I had upgraded GDAL to 2.3.1 in these wheel
builds, and was blaming my inability to access these jp2 files on a lack of
credentials, but maybe something has changed between Rasterio and JP2 in
2.3.0? I'm going to revert to 2.2.4, the version in the 1.0.2 wheels.
Sean,

Your wheels use openjpeg for GDAL JPEG2000 capability, right ? There has been
a change in GDAL ./configure in 2.3.0 related to openjpeg detection. It now
requires pkg-config to be available to detect openjpeg.

Even

--
Spatialys - Geospatial professional services
http://www.spatialys.com

Re: Fast X-array rasterization

Philipp Kats
 

oh sure!
here is the link:

And indeed, I am using geoseries.apply, but it is pretty much the same as looping - and obviously converting N 2d rasters into 3d will take some time as well - so I wonder if
there is a better solution here...




On Mon, Sep 10, 2018 at 11:43 AM Brendan Ward <brendan.charles.ward@...> wrote:
Hi Philipp,

How are you looping over each value to create a new raster?  Are you able to share a gist of what you are doing?