41 votes
Accepted

how to take cell size from raster using python or GDAL or RASTERIO

With PyQGIS ras = QgsRasterLayer("raster.tif") pixelSizeX= ras.rasterUnitsPerPixelX() pixelSizeY = ras.rasterUnitsPerPixelY() print pixelSizeX 2.11668210081 print pixelSizeY 2.11685012701 With GDAL ...
user avatar
  • 51.2k
31 votes
Accepted

Using Rasterio or GDAL to stack multiple bands without using subprocess commands

Using rasterio you could do import rasterio file_list = ['file1.tif', 'file2.tif', 'file3.tif'] # Read metadata of first file with rasterio.open(file_list[0]) as src0: meta = src0.meta # ...
user avatar
26 votes
Accepted

numpy array to GTiff using rasterio without source raster

This ended up being more straightforward than I thought, with all of the capabilities lying in the rasterio.open function. Here is an example using a proj4 string instead of wkt. import rasterio ...
user avatar
  • 1,034
21 votes
Accepted

How can I superimpose a geopandas dataframe on a raster plot?

For me it works if I pass the matplotlib ax object explicitly to rasterio.plot.show: fig, ax = plt.subplots(figsize=(15, 15)) rasterio.plot.show(raster, ax=ax) countries.plot(ax=ax, facecolor='none', ...
user avatar
  • 3,428
18 votes
Accepted

Getting pixel values at single point using rasterio

The Python API method that supports the rio-sample command is documented here: https://rasterio.readthedocs.io/en/latest/api/rasterio._io.html#rasterio._io.DatasetReaderBase.sample src.sample() takes ...
user avatar
  • 8,538
18 votes
Accepted

Finding pixel location in raster using coordinates?

A rasterio way of doing this is pretty simple. Note this requires your raster be in the same projection as your coordinates. You can of course project your coordinates on the fly, but that's another ...
user avatar
  • 52.8k
18 votes

Make shapefile from raster-bounds in Python

With shapely box and GeoPandas import rasterio as rio ra = rio.open("raster.tif") bounds = ra.bounds Convert bounds to shapely geometry from shapely.geometry import box geom = box(*bounds) print(...
user avatar
  • 51.2k
17 votes
Accepted

Reading raster files by block with rasterio?

Reading Rasters by block can be done in rasterio and I'd argue it's easier than in GDAL. There is even a tutorial on windowed read/write over at GitHub. Let's take a look at the read functions ...
user avatar
  • 9,460
17 votes
Accepted

Rasterio: Use case for src.meta vs. src.profile when reading and writing rasters?

The meta property contains the basic raster metadata. The profile property is a super set of meta which includes dataset creation options (i.e inc. tiling, block size, compression etc...). You would ...
user avatar
  • 52.8k
16 votes
Accepted

Converting NetCDF dataset array to GeoTiff using rasterio Python

I would recommend looking into rioxarray for your dataset. You can open your dataset like so: import rioxarray import xarray xds = xarray.open_dataset('D:\Weather data\et_01012018.nc') If your CRS ...
user avatar
  • 5,827
16 votes
Accepted

Creating an in memory rasterio Dataset from numpy array

You need to recreate a DatasetReader manually and you can use a MemoryFile to avoid writing to disk. You can re-use the metadata from the input raster in the DatasetReader, but you'll need to modify ...
user avatar
  • 52.8k
16 votes
Accepted

Get a window from a raster in rasterio using coordinates instead of row/column offsets, width and height

Use the rasterio.windows.from_bounds function. No need for a class or mixin. import rasterio from rasterio.windows import from_bounds from rasterio.enums import Resampling with rasterio.open(...
user avatar
  • 52.8k
15 votes

Using Rasterio or GDAL to stack multiple bands without using subprocess commands

If using GDAL 2.1+ it's as simple as gdal.BuildVRT then gdal.Translate: from osgeo import gdal outvrt = '/vsimem/stacked.vrt' #/vsimem is special in-memory virtual "directory" outtif = '/tmp/stacked....
user avatar
  • 52.8k
13 votes
Accepted

Plot shapefile on top of raster using plot and imshow from matplotlib

I added that recipe to the rasterio documentation. Since it was such a simple shape, in this case I just unzipped the coords in the single record contained by the shapefile. That is, x, y = zip(*...
user avatar
  • 812
13 votes

Using GeoPandas with rasterstats for zonal statistic gives error

For anyone else that stumbles upon this you can do: To get just one stat: gdf['mean'] = pd.DataFrame( zonal_stats( vectors=gdf['geometry'], raster='raster.tif', stats='...
user avatar
13 votes

Rasterio: tool for creating mosaic?

import rasterio import numpy as np from rasterio.merge import merge from rasterio.plot import show src1 = rasterio.open('/path/to/your/raster1') src2 = rasterio.open('/path/to/your/raster2') # ...
user avatar
  • 1,614
12 votes
Accepted

Extracting data from a raster

You can use rioxarray for this: import rioxarray from pyproj import Transformer # convert coordinate to raster projection lon = 15.174866 lat = 43.169129 rds = rioxarray.open_rasterio("C:/path-to-...
user avatar
  • 5,827
12 votes
Accepted

Displaying GeoTIFF file in Python

Yes. It is just an issue related with the parameters. Without vmin and vmax, matplotlib uses nodata value as minimum value, and stretches the color map. So, you need to specify min and max value: ...
user avatar
11 votes
Accepted

Crop a raster using rasterio and geopandas

The issue is resolved. The issue was I misread the documentation. On a second read, the rasterio.mask documentation clearly states that polygons should be a list of GeoJSON-like dicts. I found the ...
user avatar
  • 371
11 votes

Finding pixel location in raster using coordinates?

Another way to do this is to use the rasterio.transform.rowcol() method described in the rasterio transform docs. Example: import numpy as np import rasterio xs = np.array([130.5, 146.0]) ys = np....
user avatar
  • 233
11 votes

Creating an in memory rasterio Dataset from numpy array

I wrote a wrapper for rasterio.mask.mask that accepts numpy arrays as inputs. def mask_raster_with_geometry(raster, transform, shapes, **kwargs): """Wrapper for rasterio.mask.mask ...
user avatar
  • 111
10 votes
Accepted

Calculating NDVI with Rasterio

First I would use bands 4(red) and 5(nir) for Landsat 8 according to the description of the OLI instrument, and 3(red) and 4(NIR) for the Landsat TM and ETM. Second, you define an output in dtype=...
user avatar
  • 48.6k
10 votes
Accepted

Reading S3 Sentinel-2 image files with rasterio

I believe you are using the wrong bucket name. It should be s3://sentinel-s2-l1c, not s3://sentinel-pds Try: url = 's3://sentinel-s2-l1c/tiles/10/S/DG/2015/12/7/0/B01.jp2' See more info here: http:/...
user avatar
  • 754
10 votes
Accepted

Cannot install a version of GDAL via Anaconda that permits reading BigTIFFs

Issues: The Anaconda default gdal may be built without BigTIFF support. If I create a non conda-forge env, i.e conda create -n testgdal gdal I can reproduce the md['DMD_CREATIONOPTIONLIST'].find('...
user avatar
  • 52.8k
9 votes
Accepted

Excluding nodata value in band calculation with Rasterio?

You can get a Numpy masked array that covers up nodata values from Rasterio by adding a keyword argument: src.read(1, masked=True). Operations on a masked array do not use the covered up elements. If ...
user avatar
  • 8,538
9 votes
Accepted

How to split multiband image into image tiles using Rasterio?

Below is a simple example (rasterio 1.0.0 or later, won't work in 0.3.6). There might be better/simpler ways (and there is an easier way if your raster is internally tiled and the tile block sizes ...
user avatar
  • 52.8k
9 votes
Accepted

rasterio "invalid dtype: 'bool'"

Yes, you can write a one bit raster with rasterio*. You need to: write to a format that supports a 1bit dataset, such as GeoTIFF; ensure your numpy array is np.uint8/ubyte so rasterio doesnt raise ...
user avatar
  • 52.8k
9 votes
Accepted

Convert GeoJSON to GeoPandas GeoDataframe

With a list of GeoJSON-like Python geo interface geometries, simply use shapely.geometry.shape (GeoPandas uses shapely, see also Python Geo_interface applications) from shapely.geometry import shape ...
user avatar
  • 51.2k
8 votes
Accepted

rasterio masking with GeoJSON from file issue

If your GeoJSON geometry is unique (not a "type": "FeatureCollection") as in geoms = {'type': 'Polygon', 'coordinates': [[(250542.40328375285, 141691.07089614146), (250641.30366207045, 141400....
user avatar
  • 51.2k
8 votes
Accepted

Update geotiff metadata with rasterio

The previous answers are misleading or wrong. To modify the nodata value of a GeoTIFF with Rasterio, do this with rasterio.open(tiffname, 'r+') as dataset: dataset.nodata = -32767 The project ...
user avatar
  • 8,538

Only top scored, non community-wiki answers of a minimum length are eligible