Date   
Re: reading raster with masked=True leads to error

Eyal Saiet
 

This is quite embarrassing, I cannot believe I did not catch that myself.

The link is very useful thanks!


On Fri, Nov 15, 2019 at 12:32 AM Guillaume Lostis <g.lostis@...> wrote:

Hi,

I think the reason you get an error is that you have misspelled the argument (you wrote maksed instead of masked).

FYI, the list of arguments of the read() function can be found in the documentation, here: https://rasterio.readthedocs.io/en/stable/api/rasterio.io.html#rasterio.io.DatasetReader.read

Best,

Guillaume Lostis



--


Eyal Saiet


The mind is not a vessel to be filled, but a fire to be kindled. Plutarch

Re: reading raster with masked=True leads to error

Guillaume Lostis
 

Hi,

I think the reason you get an error is that you have misspelled the argument (you wrote maksed instead of masked).

FYI, the list of arguments of the read() function can be found in the documentation, here: https://rasterio.readthedocs.io/en/stable/api/rasterio.io.html#rasterio.io.DatasetReader.read

Best,

Guillaume Lostis

reading raster with masked=True leads to error

Eyal Saiet
 

Hello,
I wanted to try using read(1,masked=True) on an imported DEM, but I get the following error below:
image.png

Since I could not verify what arguments the read() method accepts, I am asking here.

thanks

--


Eyal Saiet

The mind is not a vessel to be filled, but a fire to be kindled. Plutarch

Rasterio 1.1.1

Sean Gillies
 

Hi all,

Rasterio 1.1.1 source distribution and wheels for macosx and manylinux1 are on PyPI now. Here are the changes from the change log:

Bug fixes:

- Calling a dataset's sample method with coordinates outside the extent of a
  dataset with a nodata value of None has raised a TypeError (#1822). Now, it
  gives the values we would get from a boundless read of the dataset's values.
- Use new set_proj_search_path() function to set the PROJ data search path. For
  GDAL versions before 3.0 this sets the PROJ_LIB environment variable. For
  GDAL version 3.0 this calls OSRSetPROJSearchPaths(), which overrides
  PROJ_LIB (#1823).
- Check for header.dxf file instead of pcs.csv when looking for installed GDAL
  data. The latter is gone with GDAL 3.0 but the former remains (#1823).
- RGB rasters are now properly displayed by rasterio.plot.show (#1650).

Notes:

- The wheels on PyPI include GDAL 2.4.3 with a patch that addresses the
  multithreading issue reported in #1828.

Some of you have been reporting trouble with concurrent reads of GeoTIFFs on S3 and I'm optimistic that we have a fix for at least one class of these problems. Please give the 1.1.1 wheels a try (or patch your own installations of GDAL) and let us know how it goes.

--
Sean Gillies

Re: reproject from a tiff to another tiff gives "destination band 1 appears to be read only"

Sean Gillies
 

Hi,

It's a good thing that you opened the destination (population) file in "r" mode. Otherwise, it could have been overwritten.

You must create an output array for the reprojected data and pass it as the second argument, along with dst_crs and dst_transform keyword arguments.

with rasterio.open("population") as population, rasterio.open("lights") as lights:

    reprojected = numpy.zeros(population.shape, dtype=lights.dtypes[0])

    reprojected, transform = rasterio.warp.reproject(rasterio.band(lights, 1), reprojected, dst_crs=population.crs, dst_transform=population.transform)

The reprojected array will be aligned with and have the same shape as lights.read(1).

On Mon, Nov 11, 2019 at 6:28 AM <simonm3@...> wrote:
Please can someone explain how I align a tiff file of lights with another of population to create two numpy arrays the same size and aligned. I tried:

    population = rio.open(population)
    lights = rio.open(lights)
    res = rio.warp.reproject(rio.band(lights,1), rio.band(population,1))

This gives "CPLE_IllegalArgError: Destination band 1 appears to be read-only."

I was expecting it to read the meta data (crs, height, width and affine transform) from source and destination then return something with the data of source but the shape of destination.



--
Sean Gillies

reproject from a tiff to another tiff gives "destination band 1 appears to be read only"

simonm3@...
 

Please can someone explain how I align a tiff file of lights with another of population to create two numpy arrays the same size and aligned. I tried:

    population = rio.open(population)
    lights = rio.open(lights)
    res = rio.warp.reproject(rio.band(lights,1), rio.band(population,1))

This gives "CPLE_IllegalArgError: Destination band 1 appears to be read-only."

I was expecting it to read the meta data (crs, height, width and affine transform) from source and destination then return something with the data of source but the shape of destination.

Re: DataReader.crs is empty with new environment?

Denise Draper
 

Yup, that is it: the environment variables are not set.  I may have set them manually before.
I'll have to investigate how conda is supposed to handle this (since the variables should only be set when you are in the right environment), but in the meantime it seems to work to set them manually.

-- 
  Denise Draper
  draperd@...


Re: DataReader.crs is empty with new environment?

Yann-Sebastien Tremblay-Johnston
 

Is your GDAL_DATA  and/or PROJ_LIB environment variables are set in your new windows 10 environment? Typically these are set for you when you create a new conda environment and install GDAL. Is it possible that they are set but are somehow not readable? Have you tried opening another image and reading the crs? does it happen only with this one specific image?

DataReader.crs is empty with new environment?

Denise Draper
 

Hi ---

I just got a new (Windows 10) laptop, and recreated my software environment from my old laptop.  Code that had previously worked didn't.  In particular, on opening a file, I find the crs field is empty:

import rasterio
foo = rasterio.open('path_to_file.tif')
foo.crs  # returns nothing

I assume it is some sort of version difference at fault.   However, I have tried upgrading on an old machine and downgrading on the new one, and neither helped.   So far I've only looked at the rasterio version and the gdal version.  What else should I be looking at?    I use miniconda to create the software environment if that matters.

Here's a summary of version combinations that I tried

Old Windows 10 machine
rasterio version 1.0.22   gdal version 2040100 ==> crs is CRS.from_epsg(32629)

Cloud Linux VM
rasterio version 1.0.28  gdal version 3000100 ==> crs is CRS.from_epsg(32629)
rasterio version 1.1.00  gdal version 3000100 ==> crs is CRS.from_epsg(32629)

New Windows 10 machine
rasterio version 1.1.00  gdal version 3000100 ==> crs is NULL
rasterio version 1.0.28  gdal version 3000100 ==> crs is NULL

Also: I've independently verified (with gdalinfo) that the data file does have a proper CRS, which is the same from all machines.

Re: setting BIGTIFF='IF_SAFER' not solving BigTiff problems

Luke Pinner
 

How did you install rasterio?

If you installed via Anaconda / conda using the default channel, your gdal may lack bigtiff support.

https://github.com/ContinuumIO/anaconda-issues/issues/9887

setting BIGTIFF='IF_SAFER' not solving BigTiff problems

Amine Aboufirass
 

Dear All,

I would like to point to this issue in which I brought up the BIGTIFF limitation with regards to merging rasters using rasterio:


It turns out that the suggested solution of adding BIGTIFF="IF_SAFER" to the rasterio.open arguments does not remove the problem. I still get an error:

rasterio._err.CPLE_NotSupportedError: A 140401 pixels x 36001 lines x 1 bands Int16 image would be larger than 4GB but this is the largest size a TIFF can be, and BigTIFF is unavailable.  Creation failed.

Is there a way to circumvent or workaround this? Why is BIGTIFF="IF_SAFER" not working?

Thanks,

Amine

Re: Clarification on rio bounds --sequence/--collection

Sean Gillies
 

Hi,

I agree, it is confusing. And there is a bug here. The --sequence/--collection switch was intended to allow printing a sequence of JSON texts (GeoJSON feature or bbox) or a single GeoJSON feature collection. But the latter feature is disabled by the reuse of --collection to determined the type of JSON output. 


On Thu, Oct 17, 2019 at 10:44 PM Yann-Sebastien Tremblay-Johnston <yanns.tremblay@...> wrote:
Hi,

In the usage description of `rio bounds --help` 

```
...

 --sequence / --collection  Write a LF-delimited sequence of texts containing
                             individual objects (the default) or write a
                             single JSON text containing a feature collection
                             object.
...

 --collection               Output as GeoJSON feature collection(s).
```

confuses me. The first instance of `--collection` seems to imply that if several input files are passed like `rio bounds --collection *.tif`, then the output would be a single GeoJSON document with a single FeatureCollection containing multiple Features corresponding to each tif file expanded in the shell. The second instance of `--collection` is perhaps more honest in setting expectation that a FeatureCollection object will be created for each tif file. I'm wondering if the first I described was ever meant to be intended behavior? or has there been a regression? 

Thank you so much!



--
Sean Gillies

Re: masking using another raster

Sean Gillies
 

Hi Matt,

Yes, I think that's the only builtin command that will do it. Enhancing rio-stack to more clearly provide this feature might be feasible. Note that few formats support addition of a mask and that adding a mask generally requires creation of a new file.


On Fri, Oct 11, 2019 at 3:12 PM Gregory, Matthew <matt.gregory@...> wrote:
Hi all,

Simple question, but I want to make sure I'm doing this the canonical way. 

Assume two single-band rasters with identical projections/transforms and one of the rasters has NODATA pixels defining an irregular shape and 1 otherwise.  Is the best way to mask the other raster to use rio calc? e.g.:

  rio calc "(* (read 1) (read 2))" -o out.tif input.tif mask.tif

This works, but didn't know if a different utility is better suited.

thanks, matt





--
Sean Gillies

Re: reading hdf4 files with rasterio

Sean Gillies
 

Rasterio doesn't improve the situation with subdatasets very much. You would still need to discover them via a rasterio dataset's subdatasets property. Something like

with rasterio.open("file.hdf4") as dataset:
    for name in dataset.subdatasets:
        with rasterio.open(name) as subdataset:
            print(subdataset.profile)


On Mon, Oct 21, 2019 at 2:39 AM Amine Aboufirass <amine.aboufirass@...> wrote:
Dear All,

Is there a way to read and write hdf4 datasets using rasterio? I am working with DAAC data, for example https://e4ftl01.cr.usgs.gov//DP107/MOLA/MYD13Q1.006/2013.01.09/MYD13Q1.A2013009.h09v07.006.2015254175244.hdf

I found some examples of reading the data using gdal but extracting all relevant subdatasets and writing to a new GeoTiff using gdal would be a headache, I think.

Please let me know your thoughts. Thanks in advance for your help.

Regards,

Amine



--
Sean Gillies

reading hdf4 files with rasterio

Amine Aboufirass
 

Dear All,

Is there a way to read and write hdf4 datasets using rasterio? I am working with DAAC data, for example https://e4ftl01.cr.usgs.gov//DP107/MOLA/MYD13Q1.006/2013.01.09/MYD13Q1.A2013009.h09v07.006.2015254175244.hdf

I found some examples of reading the data using gdal but extracting all relevant subdatasets and writing to a new GeoTiff using gdal would be a headache, I think.

Please let me know your thoughts. Thanks in advance for your help.

Regards,

Amine

Clarification on rio bounds --sequence/--collection

Yann-Sebastien Tremblay-Johnston
 

Hi,

In the usage description of `rio bounds --help` 

```
...

 --sequence / --collection  Write a LF-delimited sequence of texts containing
                             individual objects (the default) or write a
                             single JSON text containing a feature collection
                             object.
...

 --collection               Output as GeoJSON feature collection(s).
```

confuses me. The first instance of `--collection` seems to imply that if several input files are passed like `rio bounds --collection *.tif`, then the output would be a single GeoJSON document with a single FeatureCollection containing multiple Features corresponding to each tif file expanded in the shell. The second instance of `--collection` is perhaps more honest in setting expectation that a FeatureCollection object will be created for each tif file. I'm wondering if the first I described was ever meant to be intended behavior? or has there been a regression? 

Thank you so much!

Re: Occasional "not recognized as a supported file format" errors when reading from S3

Dion Häfner
 

Hey Sean,

I don't know if you've seen our update on the rasterio issue tracker:

https://github.com/mapbox/rasterio/issues/1686#issuecomment-541737735

As it turns out, `GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR` was the culprit. Do you know if that is related to a known bug in GDAL, or whether this is something that should be reported separately?

For now, we have disabled the flag in Terracotta, and everything seems to work fine, but we did take a performance hit from that. So in the long run, it would be great if that could be fixed.

Best,
Dion

On 15/06/2019 16.51, Sean Gillies via Groups.Io wrote:
Dion, Jonas:
There is a new report in the FIona tracker that may be related: https://github.com/Toblerity/Fiona/issues/761. There I advised turning on CPL_CURL_VERBOSE=YES to see more details about the HTTP requests that are failing. It will cause the volume of your logs to explode and it may be a challenge to find the signal in the noise, but I believe we will find some clues in the HTTP responses.
On Wed, Jun 12, 2019 at 1:55 AM Jonas <@j08lue <mailto:@j08lue>> wrote:
OK, sorry, forget what I said - the same map tile that earlier was
failing is now served fine. :/
--
Sean Gillies

Re: Flaw in rasterio design?

Sean Gillies
 

Hi all,

Since I haven't written much about Rasterio's design principles or goals beyond what is at https://rasterio.readthedocs.io/en/latest/intro.html#philosophy it is only natural that these questions come up.

Rasterio's stated design goals are to eliminate the biggest "gotchas" of the GDAL Python bindings and make GDAL configuration a little more Pythonic. The unstated design principles include minimizing the invention of new concepts and minimizing the risk of project failure. Rasterio's new concept of using a context manager to patch GDAL's runtime configuration has already used up most of our new invention budget, I believe. Whether or not this was the right choice depends on your perspective. It has certainly made my own team's code less brittle and more easy to reason about.

Not developing application-specific classes has been one way of minimizing the risk of project failure. Rasterio is going to stay useful for a long time by staying relatively small and behind the cutting edge. We have very limited resources and if we had used them to build a GIS application framework we would have been at risk (in my view) of missing the mark or being swept away by the new flavor of the month as we've seen in other domains. In my team's projects that use Rasterio we do have application classes, but they are constantly changing and I'm very glad not to have published them and committed myself to supporting their use, forever, for free.

On the other hand, count me in for helping design raster data protocols and standards. I'm all for them, but I do not think Rasterio is the right place to experiment with them.


On Tue, Oct 15, 2019 at 3:33 AM Dion Häfner <dion.haefner@...> wrote:
There are some libraries that try and implement a raster datatype, which
rasterio explicitly does not. Strangely, none of them is as successful
as rasterio.

I see rasterio as a tool to get the job done with minimal abstraction,
and I have the highest respect for this type of design. I like not
having to wrap my head around yet another container type. But others can
speak more about the original reasoning behind that design choice, which
I don't know.

I agree that in an ideal world with well-funded open source development
we would probably have pandas or xarray for raster data that doesn't
leak GDALisms into its abstractions. Until then, I suggest you come up
with a more concrete suggestion on how to improve the tools we have.

Best,
Dion

On 15/10/2019 11.14, Amine Aboufirass via Groups.Io wrote:
> Hi Dion,
>
> As an experienced user you find yourself juggling between dataset
> objects and arrays/transform objects way too often. It slows down the
> workflow and makes rasterio a pain to use. The choice of where to burn a
> file or not should be left up to the user. Dataset objects should be
> abstract and in memory at all times unless the user wishes for them to
> be written as final output.
>
> I feel like the attachment to these file abstractions is an artefact
> from whatever rasterio is built on top of (GDAL? C++?) and is not
> pythonic at all. Take for instance the pandas dataframe object. All
> operations are based on the in-memory object, or a pointer thereto. The
> choice of whether to write to excel, csv or heaven-knows-what other kind
> of file format is only a method of the dataframe object and not some
> kind of intrinsic part of it.
>
> In python the abstraction should represent the physical file and not the
> opposite.
>
> On Tue, Oct 15, 2019 at 11:05 AM Dion Häfner <dion.haefner@...
> <mailto:dion.haefner@...>> wrote:
>
>     Hey Amine,
>
>     usually there are enough functions to allow you to have your cake and
>     eat it, too.
>
>     Instead of reproject, I usually find myself using a VRT. An alternative
>     to mask that returns an ndarray is raster_geometry_mask. Some of this
>     you will have to roll yourself, though.
>
>     I think most of the functions that burn rasters directly (such as mask)
>     are meant for less experienced users that just want to get the job
>     done.
>     Maybe there should ideally always be two versions of each function, one
>     operating at dataset level, one at data (+ metadata) level. But as an
>     experienced user you can pretty much always work around that limitation.
>
>     Hope that helps,
>     Dion
>
>     On 15/10/2019 10.51, Amine Aboufirass via Groups.Io wrote:
>      > Dear All,
>      >
>      > Does anyone have any difficulty with input and output datatypes for
>      > rasterio functions in general? Many of the python functions in the
>      > library have inconsistent input and output. For instance,
>     sometimes they
>      > have ndarrays and transforms as input/output and sometimes they have
>      > physical datasets as input output. This is very annoying when
>     trying to
>      > develop a generic workflow. Here is a prime example:
>      > *
>      > *
>      > *rasterio.warp.reproject * takes ndarrays and source/destination
>      > transforms to do its work while *mask.mask *takes a dataset
>     opened in
>      > "r" mode. So then what if I first want to open a dataset, clip
>     and then
>      > reproject in my script? I would have to open the file as a
>     dataset in
>      > "r" mode, clip and then burn it into a file. Finally I would have to
>      > open this new file, extract its contents and transform and pass that
>      > into the *reproject *function and burn a third file.
>      >
>      > At first I thought memory files were a good solution, but in fact
>     they
>      > are not because we still see many functions which require
>     non-dataset
>      > formats to do their work. I really think this is a flaw in rasterio
>      > design. There should be one object type which can be passed into
>      > different functions regardless of what they do.
>      >
>      > In the meantime how can I adapt my code so I don't have to be
>     aware of
>      > input and output datatypes?
>      >
>
>
>
>





--
Sean Gillies

Re: Flaw in rasterio design?

Amine Aboufirass
 

Dear Alan, 

I just had a look. It looks like they are in fact taking those things into consideration. Very exciting. I can't wait to put this to use.

On another note, such a library would not have been possible without rasterio, so where rasterio shines (IMHO) is as a lower-level library upon which more pythonic libraries are constructed.

Regards,

Amine

On Tue, Oct 15, 2019 at 3:09 PM Amine Aboufirass via Groups.Io <amine.aboufirass=gmail.com@groups.io> wrote:
Hi Alan,

Thanks for your response. I am now a bit more intrigued by xarray, but I have a couple of questions:
  • Wouldn't you then lose all geospatial information? Things like affine and CRS are really important
  • Would I then be able to pass such an xarrayish object to rasterio functions or do similar things based on the rasterio package? I still need access to functions like masking, reprojecting, resampling and warping.
I'll take a look myself as well to check if it can reduce my suffering :P.

Regards,

Amine

On Tue, Oct 15, 2019 at 3:05 PM Alan Snow <alansnow21@...> wrote:
Hi Amine,

rasterio is definitely a great tool. But, if you are looking for an interface to it similar to pandas, you may be interested in rioxarray. It wraps the rasterio code and gives you an xarray (n-dimensional pandas) interface.

Here is a link to some examples that show how to do operations such as reproject and clip.
https://corteva.github.io/rioxarray/html/examples/examples.html

Hope this helps,
Alan

Re: Flaw in rasterio design?

Amine Aboufirass
 

Hi Alan,

Thanks for your response. I am now a bit more intrigued by xarray, but I have a couple of questions:
  • Wouldn't you then lose all geospatial information? Things like affine and CRS are really important
  • Would I then be able to pass such an xarrayish object to rasterio functions or do similar things based on the rasterio package? I still need access to functions like masking, reprojecting, resampling and warping.
I'll take a look myself as well to check if it can reduce my suffering :P.

Regards,

Amine

On Tue, Oct 15, 2019 at 3:05 PM Alan Snow <alansnow21@...> wrote:
Hi Amine,

rasterio is definitely a great tool. But, if you are looking for an interface to it similar to pandas, you may be interested in rioxarray. It wraps the rasterio code and gives you an xarray (n-dimensional pandas) interface.

Here is a link to some examples that show how to do operations such as reproject and clip.
https://corteva.github.io/rioxarray/html/examples/examples.html

Hope this helps,
Alan