[GRASS-user] Appropriate input/output range for i.atcorr WorldView2?

[Markus Neteler]

Hi Paige,

On Tue, Apr 7, 2020 at 1:45 AM Paige Byassee <byasseepaige at gmail.com> wrote:
>
> Can anyone speak to the appropriate range/datatype of WorldView2 rasters when using the i.atcorr module?
>
> I.atcorr requires the input raster be converted to at-sensor-radiance first, which I have done according to DigitalGlobe documention.
> DG suggests all transformations be done in FLOAT32 datatype, and i.atcorr requires a specified input/output range that is positive.

Yes, as per manual
https://grass.osgeo.org/grass78/manuals/i.atcorr.html:

Input range
  range=min,max
  Default: 0,255

.. the default range is 8bit. Since you use WorldView2, the range
parameter needs to be modified (see below).

I have a WorldView2 scene available here, so I am looking into the
details with gdalinfo:


Metadata:
  AREA_OR_POINT=Area
  METADATATYPE=DG
  TIFFTAG_COPYRIGHT=(C) COPYRIGHT 2018 DigitalGlobe, Inc., Longmont CO USA 80503
  TIFFTAG_DATETIME=2018:12:20 09:13:25
  TIFFTAG_IMAGEDESCRIPTION={
  bandList =
  [
    6;
...
    9;
  ]
}
  TIFFTAG_MAXSAMPLEVALUE=2047
  TIFFTAG_MINSAMPLEVALUE=0
Image Structure Metadata:
  INTERLEAVE=BAND
RPC Metadata:
...
Band 1 Block=3217x162 Type=UInt16, ColorInterp=Red
  Overviews: 1609x1660, 805x830, 403x415, 202x208
Band 2 Block=3217x162 Type=UInt16, ColorInterp=Green
...

So, UInt16 also here. Now:

> I convert from the original DN's (delivered in UINT16) to radiance (FLOAT32),

I am not entirely convinced that you need to convert UINT16 to FLOAT32
like that.
I believe that you need to apply the formula along with the respective
band calibration values ABSCALFACTOR and  EFFECTIVEBANDWIDTH (the same
is also stored in the IMD files) from the XML metadata files:

grep 'ABSCALFACTOR\|EFFECTIVEBANDWIDTH'
058891334020_01_P001_MUL/18DEC03184338-M2AS-058891334020_01_P001.XML
            <ABSCALFACTOR>9.295654000000000e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>4.730000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>9.748051000000001e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>5.430000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>7.541495000000000e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>6.300000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>5.101088000000000e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>3.740000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>1.103623000000000e-02</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>5.740000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>4.539619000000000e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>3.930000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>1.224380000000000e-02</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>9.890000000000000e-02</EFFECTIVEBANDWIDTH>
            <ABSCALFACTOR>9.042234000000000e-03</ABSCALFACTOR>
            <EFFECTIVEBANDWIDTH>9.959999999999999e-02</EFFECTIVEBANDWIDTH>

Citing from "Absolute Radiometric Calibration, Prepared By: Michele A. Kuester"
https://dg-cms-uploads-production.s3.amazonaws.com/uploads/document/file/209/ABSRADCAL_FLEET_2016v0_Rel20170606.pdf

  "The top-of-atmosphere radiance, L, in units of Wμm-1m-2sr-1, is
then found from the DigitalGlobe
    image product for each band by converting from digital numbers
(DN) using the equation,

    L = Gain * DN * (abscalfactor/effective bandwidth) + Offset

    The TDI specific "abscalfactor" and "effectiveBandwidth" are
delivered with the imagery in the metadata file.
    The digital number, DN, is the pixel value found in the imagery.
The Gain and Offset are the absolute
    radiometric calibration band dependent adjustment factors that are
given in Table 1. Note that these
    are not necessarily stagnant values and they are revisited annually.
 "

(the "Table 1" is found in the same PDF file above. You may want to
check if a newer table version exists).

You can apply this formula with r.mapcalc, for each band to obtain TOA
data from DN.

> and while my radiance values are positive, the datatype range is not. In another sense, WV2 is encoded with an 11-bit dynamic range (i.e. 2^11-bit range), but I am unsure if any of these values are right for the module... For example, Sentinel 2 L1C are also encoded as UINT16 with 15 significant bits. However, documentation states that the correct i.atcorr Sentinel range is 1-10000; 10000 being a specified QUANTIFICATION_VALUE by Sentinel.
>
> Is there a WorldView 2 equivalent to this quantification value? Is there an alternative range specific to WV2? Might this range be affected by the radiance conversion? Respectively, what would be the appropriate output/rescale range?

I hope the explanations above help.
Please let us know how it goes.

> In general- Is the i.atcorr module seeking the literal range of the image pixel values, t
> he range of the datatype they're in, or the range of data the satellite can collect?

Well, unfortunately i.atcorr isn't yet clever enough to parse the
values itself from the metadata. So, one needs to find the values and
apply them.

Here is an overview:
https://grasswiki.osgeo.org/wiki/Atmospheric_correction

and this page should IMHO be updated (your input is most welcome here
- just register there to edit or comment here)
https://grasswiki.osgeo.org/wiki/WorldView#Spectral_Radiance


Best,
Markus

--
Markus Neteler, PhD
https://www.mundialis.de - free data with free software
https://grass.osgeo.org
https://courses.neteler.org/blog