Hot answers tagged

23

I have finally gotten around to improving this function. I found that for my purposes, it was fastest to rasterize() the polygon first and use getValues() instead of extract(). The rasterizing isn't much faster than the original code for tabulating raster values in small polygons, but it shines when it came to large polygon areas that had large rasters to be ...


19

I ran a test to determine how the speed and quality differs between the two methods, here are the results: Input data 4-band NAIP DOQQ image in .img format (349.34MB) A feature class used as the mask/clipper Performance Three trials were performed and benchmarked. The Clip (Data Management) method is significantly faster than the Extract by Mask (...


15

You can use a conditional statement. The issue with previous recommendations is that when you rasterize your polygons (which is necessary) the background, that does not contain polygons, will be NoData resulting in corresponding areas in the output also being NoData. You will need to set your analysis extent to your original raster and then set a background ...


12

The documentation is indeed a little confusing. It states that: [...] all cells that are not covered by the Spatial object are set to updatevalue But actually, here covered means only if the cell's centroid is in the polygon. Indeed, mask() calls rasterize(), which states: For polygons, values are transferred if the polygon covers the center of a ...


11

Not a perfect solution but you could make use of the Geometry Generator which adds a visualised line to represent the intersection. You could then set this to overlap the original line feature. Add a new symbol layer by clicking the plus sign and select the Geometry generator as symbol layer type. Set the geoemtry type to LineString / MultiLineString and ...


10

The simple way to do this in QGIS is to use the Raster Calculator (Raster->Raster Calculator). You have a couple of options. The easiest to explain/understand is to make a unitary raster from your mask (all data set to either 1 or NoData) and then multiply your clip layer by the unitary mask layer. To ensure the extents match the mask layer, in the ...


10

Solution for single band raster. Convert polygons to raster. Use raster calculator Con(IsNull( pgonRaster), sourceRaster)


10

A raster is an array and is always square (rectangular), regardless of where the data occurs. Unless you have a perfectly square extent of data you will have NA values, this is part of the point of storing data in an array format. Now, one thing that does occur is that you have an extent much larger than necessary resulting in spurious NA values. In this ...


9

This workflow deals first with the junk floating around the edges and then solves the problem. Create a binary indicator of the "border" area, which I take to include all surrounding NoData cells. It is convenient to use any value to indicate the border and NoData for the rest, as in SetNull(Not(IsNull('X')), 1): Regiongroup the result and select the ...


8

You are looking for the gdal.RasterizeLayer function. You could then use ReadAsArray to turn the rasterized polygon into a numpy array. Based on your NetCDF file extent and rows/columns, the following code should generate you a numpy 0-1 mask that matches the NetCDF exactly. shapefile=r'whatever your shapefile path is' xmin,ymin,xmax,ymax=[139.8,-39.2,150....


8

If your GeoJSON geometry is unique (not a "type": "FeatureCollection") as in geoms = {'type': 'Polygon', 'coordinates': [[(250542.40328375285, 141691.07089614146), (250641.30366207045, 141400.7504307576), (250421.17056194422, 141512.41214821293), (250542.40328375285, 141691.07089614146)]]} and you try with rasterio.open("a.tif") as src: out_image, ...


7

Because the Raster Calculator is a spatial analyst tool, you can utilize the Mask environment. From there, you can use a variety of commands to perform the reclassification: common ones include Con, Pick, Is Null and Set Null, based on your needs. To check if a specific spatial analyst tool honors the Mask environment, simply scroll down to the bottom ...


6

If you're using Python I'd recommend using the GDAL library, which has it's own Python bindings. Assuming you've got both GDAl (see this GIS StackExchange question for details on how to install on windows) and numpy installed, your code could look something like: from osgeo import gdal import numpy as np #Open our original data as read only dataset = gdal....


6

You can create a collar using data frame or map frame within the layout view, as follows Set the data map position for both X and Y to 0 (Zero) Increase the thickness of the frame to desired value (I used 20 mm), The thickness increased from the center of the line (both outside and inside of the line), that is why I set the X and Y position to 0 to exclude ...


6

You can use Raster > Mask > Land/Sea Mask and choose Use vector as mask and the corresponding cloud band in the dropdown menu.


5

As you've probably found, removing clouds from Landsat imagery is not a trivial problem. There is a decent amount of scholarly research about methods for implementing masking and correction for clouds and their shadows. So due to the complexity of this issue, you probably won't find any one size fits all solutions that work perfectly out of the box. That ...


5

You can put the colored polygons on top, with a layer blending mode set to darken Below, have the building layer with the polygon fill in white. At the bottom, add a new layer containing one large black polygon. Without the black background:


4

Perhaps the easiest solution is to call gdal_rasterize from within your Python code (use subprocess.call or subprocess.check_call), either using your NetCDF file as the destination, or create a separate image file (GeoTIFF is always a good bet) and load it into a numpy array. It may be possible to use GDAL's In Memory Raster format, but I'm not sure what ...


4

you should gieve a try to "mask" plugin: http://plugins.qgis.org/plugins/mask/ You can either select by hand or query your objects. It takes them, dissolve them and create a hole in a square feature that is 4 times larger than select objects. It adds a semi-transparent memory layer to your project. If Memory layer saver plugin is installed, your project ...


4

Try using Clip to Shape under Data Frame properties. I use that option when I want to clip feature by specifying area of interest. It also clips the label too.


4

You can rasterize (ArcToolbox > Conversion Tools > To Raster) your polygons, and then merge the two rasters with Spatial Analyst Tools > Math > Logical > Over. While rasterizing polygons, in the Polygon To Raster dialog window, you need to use Environments -> Processing Extent Snap Raster option, to get the cells correspondent to your initial raster.


4

It works this way: Create a clip polygon for the outer border Create a clip poygon for the inner border Create a multipolygon with the second layer as a hole using Vector -> Geoprocessing -> Difference Clip the Raster to that multipolygon using Raster -> Extraction -> Clipper


4

Assuming that the polygons do not overlap in any way, you can use Tabulate Areas tool in Spatial Analyst to calculate the areas in batch. Your polygons would be the ZoneRas, and your reclassified raster would be the ClassRas: Note that if the zone layer is a polygon feature layer, this tool rasterizes it before processing, which is why it asks for a ...


4

So... my solution was as follows. Thank you Michael for the response and to others who might have other solutions. import pyproj from shapely.geometry import Polygon, Point # WGS84 datum wgs84 = pyproj.Proj(init='EPSG:4326') # Albers Equal Area Conic (aea)...


4

Speed up extracting raster (raster stack) from point, XY or Polygon Great answer Luke. You must be a R wizard! Here is a very minor tweak to simplify your code (may improve performance slightly in some cases). You can avoid some operations by using cellFromPolygon (or cellFromXY for points) and then clip and getValues. Extract polygon or points data from ...


4

I would recommend using Minimum Bounding Geometry (Data Management) followed by Clip (Data Management) with the clipping geometry option or Extract by Mask (Spatial Analyst).


4

If you mask any raster you probably will obtain NA values, even more using SpatialPolygonDataFrame as mask. You have two posible options, if SpatialPolygonDataFrame is a rectangle, use crop() before mask to reduce raster's extend. Second option, change NA values to other value, such 0 or -9999: r <- raster(ncol=10, nrow=10) m <- raster(ncol=10, nrow=...


4

You can try crop and mask instead of only masking over the Raster Layers. Try the reproducible and commented code below. In this example, using crop reduced ~ 92% the size of the objects in the R environment (memory usage). To measure the memory usage I ran the function in this post: Tricks to manage the available memory in an R session. And if you want to ...


4

Here's a version with a few tweaks: library(sp) library(raster) library(rgdal) library(rmapshaper) library(tidyverse) nx <- 361 ; ny <- 181 xmin <- -30.0 ; ymin <- 25.0 dx <- 0.25 ; dy <- 0.25 lat <- seq(ymin, ymin+(ny-1)*dy, dy) lon <- seq(xmin, xmin+(nx-1)*dx, dx) # make an empty grid instead so NA = Ocean m <- matrix(NA, ...


4

I am assuming you have access to the spatial analyst extension license. You can rasterize your polygons with Polygon to Raster, use OID/FID as value field but it doesn't matter what field you use as long as it is numeric, but before you do set your environments cell size, output extent and snap raster to your existing NDVI raster. Your polygon raster will ...


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