Extracting UINT16 from MODIS

[conradb]

Hello,
I have been using polar2grid to extract and reproject MODIS data.However, the data produced is always of type Byte. Is there a function which is able to extract the MODIS data as original outputs in UINT16?

Thanks in advance for your help.
Have a great day,
Conrad

[kathys]

Hello Conrad,

Can you provide us with the exact Polar2Grid command that you are using for this extraction? We can then advise you on how to get the output in the format that you want.

Kathy

[conradb]

Hi Kathy,

here is an example that I tested:
polar2grid.sh -r modis_l1b -w geotiff --dtype uint16 --no-enhance --fill-value 0 --vis-products --grid-configs ...

I have also tried using the command option to extract specific bands only but I get the same result.

Thanks,
Conrad

[davidh]

Hi Conrad,

Just to be sure, what version of Polar2Grid are you running?

Could you send along the log file that P2G produced from this uint16 execution? If possible, a copy of the file it produced may also speed up my debugging. Thanks.

Dave

[conradb]

Hi Dave,
I'm using v3.0. Since you asked you for the logs, I also started looking into them. I am sending you two logs. The first was based on my own grid definition (which I also attach) and works when producing 'true-color' outputs. The second, I used the default polar2grid projection.

The first bands seems to be empty but the second has something (but not everything).
Please note that the system did not allow me to upload the .log or .yaml files so I added the .tif extension to the text files. They should be readable.

Thanks and have a great day,
Conrad

[davidh]

Hi Conrad,

I see two attachments, but both seem to be valid geotiffs. Both of them also seem to have UInt16 data types. Looking at the vis01 image in python, I see it is indeed all 65535 for the swath portion of the image (0 otherwise). I'm hoping the log can tell me a bit more.

Dave

[conradb]

Hi Dave,
sorry for the silence. The output is indeed UINT16 but I'm also trying (now) to understand why there are differences in the outputs and which solution I should choose. There are a number of output options but they all seem to differ in the values produced.

1. polar2grid.sh -r modis_l1b -w geotiff --dtype uint16 --no-enhance --vis-products ->
Band 1 Block=256x256 Type=UInt16, ColorInterp=Gray
STATISTICS_MAXIMUM=65535
STATISTICS_MEAN=39055.995523145
STATISTICS_MINIMUM=0
STATISTICS_STDDEV=31095.578385324
=> seems correct but

2. polar2grid.sh -r modis_l1b -w geotiff --no-enhance --vis-products ->
Band 1 Block=256x256 Type=Float32, ColorInterp=Gray
STATISTICS_MAXIMUM=25.901718139648
STATISTICS_MEAN=4.9152803946267
STATISTICS_MINIMUM=0.0088238604366779
STATISTICS_STDDEV=4.4809508808444
=> I would have expected the same stats even though the output type is float. Why the difference?

3. polar2grid.sh -r modis_l1b -w geotiff --method nearest --vis-products ->
Band 1 Block=256x256 Type=Byte, ColorInterp=Gray
STATISTICS_MAXIMUM=175
STATISTICS_MEAN=32.773020568035
STATISTICS_MINIMUM=0
STATISTICS_STDDEV=33.222786151543
=> expected because default output is byte. But why did test 2 above produce float? What does the -no-enhance flag actually do?

4. polar2grid.sh -r modis_l1b -w geotiff --no-enhance --method nearest --vis-products ->
Band 1 Block=256x256 Type=Float32, ColorInterp=Gray
STATISTICS_MAXIMUM=26.473184585571
STATISTICS_MEAN=4.9154697913181
STATISTICS_MINIMUM=0
STATISTICS_STDDEV=4.4888636965794
=> similar to test 3 - expected as the only difference will be in the interpolation between nearest neighbour and ewa

So which is the best to keep as close to the original 12-bit values of the MODIS data?

Thanks for your help,
Conrad

[davidh]

When you specify `--no-enhance` and don't force `--dtype` the geotiff writer will try to preserve the data type of the data being saved. So for these visible bands that is 32-bit floats.

For the unexpected 16-bit values, the geotiff writer is scaling the input floating point data (the visible reflectances) to the output data type of the file (16-bit unsigned integers) and assumes that the data is already normalized to a 0-1 range. If the data had been integer data (ex. cloud type) then it would not have done this. So the reflectances which are represented in Satpy/Polar2Grid percentages from about 0 to about 100 are being multiplied by 65535 and then clipped to the valid range of the 16-bit integers (0 and 65535). Since the reflectances seem to have a maximum 25.9% and a mean of ~4% you're getting mostly clipped values (4 * 65535 is clipped to 65535). I think I'm most surprised by the low maximum value in the 32-bit float geotiff only be 25.9%. Are there no clouds in the image?

So your options are:

1. Use the 32-bit float geotiff with --no-enhance.
2. Add custom "enhancements" to scale the visible data in a predictable linear way (default is square root) and save it to --dtype uint16.

There is no way currently to get the uncalibrated values stored in the MODIS L1b files using Polar2Grid if that is what you were after. These can be loaded with Satpy by telling it to load the calibration "counts". Depending on exactly what you need this may be the best option to fully customize what you get. I'm not entirely sure what you need but I personally would choose the 32-bit float geotiff. If I was doing further processing in Python I would use Satpy directly.

Let me know if I'm misunderstanding something about what you're doing.

Dave

[conradb]

Thanks David,
that does help. My final choice for running the polar2grid is the following:
polar2grid.sh -r modis_l1b -w geotiff --dtype uint16 --vis-products

and seems to output the wanted counts. I'll keep you posted if I run into any other issues.

A second question that I have is with respect to modis_l2 functionality. While it seems to be a reader that is available through polar2grid, I'm not sure what I can use it for. I was hoping to use it for extracting cloud mask data from the following HDF type: MOD35_L2.A2023324.0935.061.2023324104948.NRT.hdf

which is made available through LANCE. However, I get error that it cannot be read. Any solutions you can propose?

Have a great day,
Conrad

[davidh]

Regarding that level 1b command:

The data saved as 16-bit unsigned integers in that are the scaled version of the square root of the reflectance percentages. You would need to provide a custom enhancement (requires modifying YAML files) to scale the data without square root and in a way that the data is normalized to 0-1 so when it is saved to the 16-bit geotiff it doesn't get clipped for being above 65535.

For level 2: support for MODIS L2 is rather limited due to lack of user interest over the years. That said, MODIS L2 in Satpy should still work. Your filename however does not match the expected filenaming schemes that Satpy knows about. If you rename the file to match this pattern:

MOD35_L2.AYYYYJJJ.000.YYYYJJJHHMMSS.hdf'

it should work. In the above the first YYYYJJJ should be the observation date, the second date is the production date/time, and the 000 is the "collection" number so 000 should be fine as it isn't used for anything in Satpy but does need to be 3 digits. The long term fix would be to add this filename scheme to Satpy. If that's something you're interested in being supported then I would recommend filing a bug report with the Satpy project at https://github.com/pytroll/satpy/issues.