New answers tagged

2

Here's a very naive way using rasterio.shutil.copy to generate VRT XML for each source raster, combine them and output a stacked VRT (equivalent to gdalbuildvrt -separate). By naive, I mean it does not check any of the properties and assumes rasters are all the same shape, extent, pixel size, data type, projection etc. import xml.etree.ElementTree as ET ...


0

If you're looking for a QGIS solution, you can check the r.to.vect GRASS function under the QGIS Processing Toolbox. These links below might help: Converting raster to vector by generating center lines? https://www.youtube.com/watch?v=TQ1sD4NSpHo


5

You're writing a floating point array as an unsigned 8 bit integer. Try src = rasterio.open('raster2.tif') profile = src.profile with rio.open('test.tif', 'w', **profile) as dst: dst.write(new_raster, 1) # No need to cast to int, as profile is already float32 or if that doesn't work, enforce the casting using .astype(rio.float32) instead of your ...


5

If you're ultimately working toward a polygon (circle) intersection with your raster, there are built-in methods for calculating a numpy mask, given a polygon. You can use this question to see how Shapely can be used to get a circle of defined radius from a centre point (you will need to convert to a projected coordinate system). Then, use rasterio.features....


0

So, I got it working using rasterio. The problem was that the Affine pivot parameter is actually in pixel coordinates, and I was passing real world coordinates. So few lines of code later, this is rotating as exected! #!/usr/bin/python from optparse import OptionParser import rasterio from affine import Affine # For easly manipulation of affine matrix from ...


3

Here's a non geocube way using rasterio and fiona (could also use the higher level geopandas which uses fiona under the hood). This loops through small chunks (aka a rasterio.windows.Window) and the features in your polygon dataset, rasterizes the features just in the small chunk area, sums the chunks, writes each summed chunk to the appropriate part of ...


0

# Look at this simple code: import os import gdal path = r'inputpath_directory' # for Loop bands = [] for i in os.listdir(path): bands.append(i) print(bands) # a list of bands # To save a raster using gdal: gdal_array.SaveArray(bands, r'outputpath_name', "GTiff", ...


0

# Look At this code: import rasterio import numpy as np import matplotlib.pyplot as plt # open the bands: B4 = rasterio.open(inputpath_name) B5 = rasterio.open(inputpath_name) # First calculate your NDVI with this function: def ndvi(red, nir): red = B4.read(1) nir = B5.read(1) red = red.astype('float64') nir = nir.astype('float64') ...


1

This should work with what you require. Read about rasterio write here d = ['Blue', 'Red', 'Green'] test = 'test.tif' folder = Path('Path/to/My/Folder') with rasterio.open("test.tif") as src: for f in folder.glob("*.tif"): if any(color in f.name for color in d): with rasterio.open(f) as med: new_raster = ...


0

I would recommend using @radouxju approach to calculating NDVI in this answer. Here is an untested example based on the link you provided and radouxju's NDVI approach: import os import rasterio import numpy as np # Filepath dataset = r'C:\HY-DATA\HENTENKA\CSC\Data\Helsinki_masked_p188r018_7t20020529_z34__LV-FIN.tif' # Outfile path outpath = r'C:\path\to\...


1

rioxarray can simplify that: https://corteva.github.io/rioxarray/stable/examples/reproject_match.html https://www.earthdatascience.org/courses/use-data-open-source-python/intro-raster-data-python/raster-data-processing/subtract-rasters-in-python/ https://carpentries-incubator.github.io/geospatial-python/07-raster-calculations/index.html


1

The rasterio documentation for the mask function is not very clear here. The shape parameter must be an iterable of geometries, not a simple geometrie. The documentation of the shape parameter (i.e. "The values must be a GeoJSON-like dict or an object that implements the Python geo interface protocol...") describes the nature of the values ...


1

You can iterate all over array, if that helps? def reclassify(arr, rows, cols): for r in range(rows): for c in range(cols): temp_value = arr[r][c] if temp_value < 20: temp_value = 4 elif temp_value >=80: temp_value = 1 elif (temp_value >=50) and (...


0

from osgeo import gdal gdal.Warp(destNameOrDestDS = outputpath, # directory output srcDSOrSrcDSTab = intputpath, # directoy intput cutlineDSName = shapefile or geojson, # vector file cropToCutline = True, # Select True copyMetadata = True, # optional dstNodata = 0) # if you have values ...


Top 50 recent answers are included