1

I am masking a raster file using rasterio. I am doing this using the mask function. The result of the masking changes the value of the pixels when loaded to GIS.

import fiona
import rasterio
import rasterio.mask

with fiona.open("test_polygon.json", "r") as shapefile:
    shapes = [feature["geometry"] for feature in shapefile]

with rasterio.open("exampe.tif") as src:
    out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
    out_meta = src.meta
    
out_meta.update({"driver": "GTiff",
                 "dtype": rasterio.float32,
                 "height": out_image.shape[1],
                 "width": out_image.shape[2],
                 "transform": out_transform})

with rasterio.open("swr_masked.tif", "w", **out_meta) as dest:
    dest.write(out_image)

The example files are located HERE.

When loaded to a GIS application like QGIS, the value of the masked data changes. It is like it was changed to another data type like int16 or something. I specified the dtype in writing to be the same with the input float32. But still have the same result.

Result in QGIS: enter image description here

Improve this question
4
  • 1
    Check with gdalinfo if the original data happens to be scaled.
    – user30184
    Jan 29 at 14:31
  • What does it mean to be scaled in this manner? In your answer below, the Scale and scale_factor are the same.
    – Nikko
    Jan 30 at 8:05
  • 1
    Scale and scalefactor mean that if the raster pixel has a value 20 it must be multiplied by 0.05 for getting the intensity of radiation as watts per square meter that seems to be what the raster data is presenting. Scale is a standard TIFF metadata and scale_factor is stored as NETCDF metadata into another place of the file. In the conversion from NetCDF into TIFF the scale_factor has been converted into Scale that is the corresponding tiff tag but the original metadata scale_factor is also kept.
    – user30184
    Jan 30 at 8:32
  • Thank you. My previous question would consider the scale and offset: gis.stackexchange.com/questions/392438/… in the output file.
    – Nikko
    Jan 30 at 9:06

1 Answer 1

Reset to default
1

The original image is scaled:

gdalinfo example.tif
...
Band 1 Block=192x21 Type=Int16, ColorInterp=Gray
  NoData Value=-32768
  Unit Type: W/m^2
  Offset: 0,   Scale:0.0500000007450581
  Metadata:
    add_offset=0
    long_name=Shortwave radiation
    missing_value=-32768
    NETCDF_VARNAME=SWR
    scale_factor=0.050000001
    units=W/m^2
    valid_max=26000
    valid_min=0

QGIS knows how to apply scaling and it divides the raw pixel values by 20 but your rasterio code does not preserve the metadata about scaling.

I do not know how to handle the scaling with rasterio but GDAL seems to do it automatically

gdalwarp -ot float32 -cutline test_polygon.json -crop_to_cutline example.tif out1.tif

enter image description here

gdalinfo out1.tif
...
Band 1 Block=2x1 Type=Float32, ColorInterp=Gray
  NoData Value=-32768
  Unit Type: W/m^2
  Offset: 0,   Scale:0.0500000007450581
  Metadata:
    add_offset=0
    long_name=Shortwave radiation
    missing_value=-32768
    NETCDF_VARNAME=SWR
    scale_factor=0.050000001
    units=W/m^2
    valid_max=26000
    valid_min=0
Improve this answer

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.