PolarGeospatialCenter / imagery_utils

Other
34 stars 10 forks source link

Fix: Assign NoData to a value outside of the valid data range for outputs of pgc_ortho.py and pgc_pansharpen.py #74

Closed power720 closed 6 months ago

power720 commented 8 months ago

I have tried two (cumulative) code changes to attempt to address the NoData issue (#73) for the subset of test images at V:\pgc\data\common\quickbase\2204_UW_Andresen_Alaska_2023\imagery_check.

Attempt 1 - Modify the WarpImage function

I modified the lib.ortho_functions.WarpImage function to assign a NoData value to the destination file based on the specified data type of the destination file.

I ran a test batch using the following command.

python pgc_ortho.py \
    --epsg=utm \
    --stretch rf \
    --outtype UInt16 \
    --dem /mnt/pgc/data/elev/dem/copernicus-dem-30m/mosaic/global/cop30_tiles_global_wgs84-height_windows.vrt \
    --format GTiff \
    --save-temps \
    /mnt/pgc/data/common/quickbase/2204_UW_Andresen_Alaska_2023/imagery_check/renamed \
    /mnt/pgc/data/common/quickbase/2204_UW_Andresen_Alaska_2023/imagery_check/ortho-16-rf-noDataTest

Note that I specified the wrong --dem VRT for the platform I was running on (WSL). The process complained that it couldn't find the reference rasters and thus the valid cell values aren't what they should be, but the NoData values behaved as expected for the output of this function.

The change successfully set the NoData value to 65535 for the intermediate images that this function produces, ortho-16-rf-noDataTest/*_warp.tif.save. However, the final images, ortho-16-rf-noDataTest/*.tif still end up with a NoData value of 0.

Attempt 2 - Modify the calcStats function

I modified lib.ortho_functions.calcStats, the function responsible for applying the stretch corrections. This is what is called on the warped image to produce the final output. I changed the hardcoded value passed to vds.GetRasterBand(band).SetNoDataValue(0) from 0 to 65535 as that is what the incoming and outgoing NoData value should be for this specific set of images.

I ran another test batch with the following command; fixing the --dem option.

python pgc_ortho.py \
    --epsg=utm \
    --stretch rf \
    --outtype UInt16 \
    --dem /mnt/pgc/data/elev/dem/copernicus-dem-30m/mosaic/global/cop30_tiles_global_wgs84-height_nunatak.vrt \
    --format GTiff \
    --save-temps \
    /mnt/pgc/data/common/quickbase/2204_UW_Andresen_Alaska_2023/imagery_check/renamed \
    /mnt/pgc/data/common/quickbase/2204_UW_Andresen_Alaska_2023/imagery_check/ortho-16-rf-noDataTest-2 

The final images have the NoData value set as 65535 as intended. However, the value of the NoData cells is modified by the stretching process, resulting in each band having a unique cell value in the NoData areas that is not 65535.

power720 commented 8 months ago

Update - NoData assignment successful

The latest commit implements the dynamic assignment of the NoData value of pansharped/ortho'd rasters based on the specified output data type (Byte, UInt16, Float32).

The change resolves the previous issue where the stretch process would modify NoData cell values. The key to mitigating this behavior was including a <NODATA> tag within the <ComplexSource> tag of the VRT. Below is an example of a VRT with this tag included.

  <VRTRasterBand dataType="Float32" band="4" blockXSize="256" blockYSize="256">
    <Metadata>
      <MDI key="NITF_ISUBCAT">816.5</MDI>
    </Metadata>
    <NoDataValue>65535</NoDataValue> <!-- Included previously. Doesn't impact stretch modifications -->
    <ComplexSource>
      <SourceFilename relativeToVRT="0">/mnt/pgc/data/common/quickbase/2204_UW_Andresen_Alaska_2023/imagery_check/ortho-16-rf-noDataTest-3/QB02_20060815222516_101001000524B400_06AUG15222516-M1BS-052800787060_01_P001_u16rf32606_warp.tif</SourceFilename>
      <SourceBand>4</SourceBand>
      <SrcRect xOff="0" yOff="0" xSize="7478" ySize="8901" />
      <DstRect xOff="0" yOff="0" xSize="7478" ySize="8901" />
      <NODATA>65535</NODATA> <!-- New addition. Tells the stretch operation to ignore NoData values -->
      <ScaleOffset>0</ScaleOffset>
      <ScaleRatio>1</ScaleRatio>
      <LUT>0:-49.8567,2047:2932.87</LUT>
    </ComplexSource>
  </VRTRasterBand>

Next steps

Finalize output data type to NoData value mapping.

power720 commented 6 months ago

Output data type to NoData value: