47
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
...
- 54.2k
36
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
...
- 1,252
32
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
# ...
- 2,262
26
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', ...
- 3,828
25
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(...
- 54.2k
23
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
with rasterio.open(filepath) as src:
rst = src.read(1, window=...
- 64k
21
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 ...
- 64k
20
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 ...
- 8,978
20
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 ...
- 7,141
19
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 ...
- 64k
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 ...
- 9,775
17
votes
Accepted
Creating an in memory rasterio Dataset from numpy array
Note: You can use rasterio.features.geometry_mask to mask your numpy array without writing a dataset (example).
Otherwise if you want to use rasterio.mask.mask, you can create a DatasetReader manually ...
- 64k
16
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....
- 64k
15
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-...
- 7,141
14
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(*...
- 894
14
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='...
- 141
14
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 ...
- 241
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')
# ...
- 1,874
13
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 ...
- 64k
13
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
...
- 54.2k
13
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:
...
- 74.7k
12
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....
- 253
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 ...
- 371
11
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 ...
- 64k
11
votes
Accepted
Merging rasters with rasterio in blocks to avoid MemoryError?
I finally solved the problem by using numpy.memmap to
create a memory-map to an array stored in a binary file on disk
and then processing the input rasters in windows and blocks. It might be slower ...
- 3,747
11
votes
Mask xarray dataset using a shapefile
Masking NetCDF time series data from shapefile using Python
You can use rioxarray. Here is an example: https://corteva.github.io/rioxarray/stable/examples/clip_geom.html
import rioxarray
import ...
- 7,141
11
votes
Accepted
Using rasterio to crop image using pixel coordinates instead of geographic coordinates
You can use a window -rasterio.windows.Window to read by pixel offsets. The georeferencing can be easily calculated from the window using the source dataset window_transform method.
import random
...
- 64k
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:/...
- 764
10
votes
Add bands' name and description to the Metadata when stacking using rasterio
If using rasterio >= 1.0, use the dataset.set_band_description(self, bidx, value) method and dataset.descriptions property.
Sets the description of a dataset band.
Parameters
----------
bidx : int
...
- 64k
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('...
- 64k
Only top scored, non community-wiki answers of a minimum length are eligible