[GRASS-user] help to use r.mapcalc with decimal/binary (modisdata)

GRASS GRASS user
1

Thank you very much Glynn and the Grass Community

This has worked perfectly!!!

Andrew

Glynn Clements glynn@gclements.plus.com
05/28/2007 02:52 AM CET

To andrew.haywood@poyry.com
cc grassuser@grass.itc.it, grassuser-bounces@grass.itc.it
bcc
Subject Re: [GRASS-user] help to use r.mapcalc with decimal/binary (modisdata)

andrew.haywood@poyry.com wrote:

Okay - unfortunately i am now totally confused. I have tried to read up on
binary data (to little or no avail!!!). My understanding is the MODIS
suf_refl_qc data that i have import is a unsigned 16 bit integer number. I
would like to select out pixels using this raster to create a cloud mask
based on some of the bit codes.

for example identify " band 1 quality = missing input" would have the
following

        • 1 0 1 1 * * * * * * * * * - where * represents placeholder

If i try and follow Glynn’s logic he is saying divide the number by 16 and
then take the modulus 16 and then equate it to “11” to get the subset.

I am unsure how I choose the “11” and how do i change the syntax to pull
out the different masks i may be interested in.

0 0 * * * * * * * * * * * * * * ideal quality all bands
0 1 * * * * * * * * * * * * * * less than ideal quality

    • 0 0 * * * * * * * * * * * * clear
    • 0 1 * * * * * * * * * * * * cloudy
    • 1 0 * * * * * * * * * * * * mixed
        • 0 0 0 0* * * * * * * * highest quality band 1
                  • 0 0 0 0* * * highest quality band 2
                    etc…

Any help would be greatly appreciated. At this stage im plugging in the
decimal numbers into binary calculator and using r.reclass to create a
MASK.

First, binary numbers are normally written with the bit 0 (the “units”
digit) on the right, so “band 1 quality = missing input” would be:

1 1 1 1 1 1
5 4 3 2 1 0 9 8 7 6 5 4 3 2 1 0

                • 1 0 1 1 * * * *

To shift everything right by N bits, divide by 2^N, so dividing by 16
(2^4) shifts everything right by 4 places, moving bits 4-7 to bits
0-3:

                        • 1 0 1 1

[Similarly, dividing a decimal number by 10^4 = 10,000 shifts it right
by 4 digits: 12345678 / 10000 = 1234]

To keep the bottom N bits and discard the rest, take the remainder
(modulus) after division by 2^N, so modulus 16 (2^4) keeps the bottom
4 bits:

0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1

[Similarly dividing a decimal number by 10^4 = 10,000 and taking the
remainder takes the bottom 4 digits: 12345678/10000 = 1234 remainder
5678, so 12345678 % 10000 = 5678]

The binary number “1011”
= 12^3 + 02"2 + 12^1 + 12^0
= 18 + 04 + 12 + 11
= 8 + 2 + 1
= 11

So, given your original bit definitions

bitsdescriptiondividemodulusresult range
0-1MODLAND QA140-3
2-3cloud state440-3
4-7band 1 data quality16160-15
8-11band 2 data quality256160-15
12atmospheric correction409620-1
13adjacency correction819220-1
14different orbit1638420-1
15unused3276820-1

For the right-hand side:

MODLAND QA (2 bits):

binarydecimaldescription
000corrected product produced at ideal quality – all bands
011corrected product produced, less than ideal quality – some or all bands
102corrected product not produced due to cloud effects – all bands
113corrected product not produced for other reasons – some or all

Cloud state (2 bits):

binarydecimaldescription
000clear
011cloudy
102mixed
113not set, assumed clear

band quality (4 bits):

binarydecimaldescription
00000highest quality
00011
00102
00113
01004
01015
01106
01117
10008dead detector; data has been copied from adjacent detector
10019solar zenith >= 86 degrees
101010solar zenith >= 85 and < 86 degrees
101111missing input
110012internal constant used in place of climatological data for at least one atmospheric constant
110113correction out of bounds, pixel constrained to extreme allowable value
111014L1B data faulty
111115not processed due to deep ocean or clouds

A shell script which uses r.mapcalc to split a QC map given as its
first argument into separate maps for each field:

#!/bin/sh
src=$1
r.mapcalc <<EOF
$src.modland = ($src / 1 ) % 4
$src.cloud = ($src / 4 ) % 4
$src.band1qual = ($src / 16 ) % 16
$src.band2qual = ($src / 256 ) % 16
$src.atmos_corr = ($src / 4096 ) % 2
$src.adj_corr = ($src / 8192 ) % 2
$src.diff_orbit = ($src / 16384) % 2
$src.unused = ($src / 32768) % 2
EOF

–
Glynn Clements glynn@gclements.plus.com

This email and any files transmitted with it are confidential, may be legally privileged and are intended solely for the use of the individual or entity to whom they are addressed. If you are not the intended recipient, you are hereby notified that any use, distribution, or reproduction of the contents of this email is strictly prohibited and may be unlawful. If you are not the intended recipient, please notify the sender by return email and destroy all copies of the original message including any attachments thereto.
Thank you.