Page 1 of 1

mapgen script with MODIS L2

Posted: Mon Jan 30, 2023 3:54 pm America/New_York
by fbecker23
Hello,

I have been using the gpt command to process MODIS and VIIRS L2 images of chl and sst. I am using SeaDAS 8.3, with OCSSW installed, and Ubuntu 22.1. In the .xml (attached at the end of the post) in summary what I do is:

read(chl image) --> BandMaths(to create chl and flags bands) -->
Collocate(to bring all the bands together) --> Reproject --> Subset --> Write(as netcdf-CF)
read(sst image) --> BandMaths(to create chl and flags bands) -->

In another post, gnwiii recommended me that the correct treatment for L2 images is to use the mapgen function so as not to generate errors on the edges of the image. I have been reading about the mapgen function (https://seadas.gsfc.nasa.gov/help-8.3.0/processors/ProcessMapgen.html), but I was not able to execute it. When I executed the command: (I'm using Python 3.7)

python path/SeaDAS/ocssw/bin/mapgen --ifile=path/A2020003181000.L2_LAC_OC.x.nc --ofile=path/mapgenTest_A2020003181000.L2_LAC_OC.x.nc --product=chlor_a,l2_flags --resolution=1km --oformat=netcdf4 --projection=mercator --verbose

I get the error:

###### l2bin ######

Running: /home/fernando/SeaDAS/ocssw/bin/l2bin flaguse=ATMFAIL,CLDICE,BOWTIEDEL l3bprod=chlor_a,l2_flags ifile=/home/fernando/Escritorio/Pruebas_gpt/chl/A2020003181000.L2_LAC_OC.x.nc ofile=mapgen_tmp.l2bin.nc area_weighting=1 prodtype=regional resolution=H
Exiting: /home/fernando/SeaDAS/ocssw/bin/l2bin flaguse=ATMFAIL,CLDICE,BOWTIEDEL l3bprod=chlor_a,l2_flags ifile=/home/fernando/Escritorio/Pruebas_gpt/chl/A2020003181000.L2_LAC_OC.x.nc ofile=mapgen_tmp.l2bin.nc area_weighting=1 prodtype=regional resolution=H returned a status 1


I need output files in netcdf format (better if netcdf-CF, that's why I'm using gpt with Snap), what kind of netcdf is the output file, CF or BEAM-DIMAP?

My goal is to get the L2 image with a regular grid (no matter the projection), the bandmath steps like the subset can be skipped. I'm looking for the regular grid so I can work the images comfortably in python. Is the procedure I'm doing to process the images using the .xml correct?

Do you know any other way to get a regular grid from a L2 image?

thank you very much for your help

------------------------------------------------------

*file.xml*

<graph id="firstGraph">
<version>1.0</version>

<node id="Read_chl">
<operator>Read</operator>
<sources/>
<parameters class="com.bc.ceres.binding.dom.XppDomElement">
<file>${filechl}</file>
</parameters>
</node>

<node id="BandMaths_chl">
<operator>BandMaths</operator>
<sources>
<sourceProduct refid="Read_chl"/>
</sources>
<parameters>
<targetBands>
<targetBand>
<name>cldice_chl</name>
<expression>l2_flags.CLDICE</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>chlfail_chl</name>
<expression>l2_flags.CHLFAIL</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>land_chl</name>
<expression>l2_flags.LAND</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>bowtiedel_chl</name>
<expression>l2_flags.BOWTIEDEL</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>straylight_chl</name>
<expression>l2_flags.STRAYLIGHT</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>orig_chl</name>
<expression>chlor_a</expression>
<description></description>
<type>float32</type>
<validExpression></validExpression>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>chl</name>
<expression>if chlor_a &gt; 8 or chlor_a &lt; 0.06 then NaN else chlor_a</expression>
<description></description>
<type>float32</type>
<validExpression></validExpression>
<noDataValue>NaN</noDataValue>
</targetBand>
</targetBands>
</parameters>
</node>

<node id="Read_sst">
<operator>Read</operator>
<sources/>
<parameters class="com.bc.ceres.binding.dom.XppDomElement">
<file>${filesst}</file>
</parameters>
</node>

<node id="BandMaths_sst">
<operator>BandMaths</operator>
<sources>
<sourceProduct refid="Read_sst"/>
</sources>
<parameters>
<targetBands>
<targetBand>
<name>cldice_sst</name>
<expression>l2_flags.CLDICE</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>bad_sst</name>
<expression>qual_sst</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>land_sst</name>
<expression>l2_flags.LAND</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>bowtiedel_sst</name>
<expression>l2_flags.BOWTIEDEL</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>straylight_sst</name>
<expression>l2_flags.STRAYLIGHT</expression>
<description></description>
<type>int32</type>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>orig_sst</name>
<expression>sst</expression>
<description></description>
<type>float32</type>
<validExpression></validExpression>
<noDataValue>NaN</noDataValue>
</targetBand>
<targetBand>
<name>sst</name>
<expression>if sst &gt; 30 or sst &lt; 0 then NaN else sst</expression>
<description></description>
<type>float32</type>
<validExpression></validExpression>
<noDataValue>NaN</noDataValue>
</targetBand>
</targetBands>
</parameters>
</node>

<node id="CollocateNode">
<operator>Collocate</operator>
<sources>
<sourceProduct refid="BandMaths_chl"/>>
<sourceProduct.1 refid="BandMaths_sst"/>>
</sources>
<parameters class="com.bc.ceres.binding.dom.XppDomElement">
<targetProductType>COLLOCATED</targetProductType>
<resamplingType>NEAREST_NEIGHBOUR</resamplingType>
</parameters>
</node>

<node id="Reproject">
<operator>Reproject</operator>
<sources>
<sourceProduct refid="CollocateNode"/>
</sources>
<parameters class="com.bc.ceres.binding.dom.XppDomElement">
<wktFile/>
<crs>PROJCS[&quot;Mercator_1SP / World Geodetic System 1984&quot;,
GEOGCS[&quot;World Geodetic System 1984&quot;,
DATUM[&quot;World Geodetic System 1984&quot;,
SPHEROID[&quot;WGS 84&quot;, 6378137.0, 298.257223563, AUTHORITY[&quot;EPSG&quot;,&quot;7030&quot;]],
AUTHORITY[&quot;EPSG&quot;,&quot;6326&quot;]],
PRIMEM[&quot;Greenwich&quot;, 0.0, AUTHORITY[&quot;EPSG&quot;,&quot;8901&quot;]],
UNIT[&quot;degree&quot;, 0.017453292519943295],
AXIS[&quot;Geodetic longitude&quot;, EAST],
AXIS[&quot;Geodetic latitude&quot;, NORTH]],
PROJECTION[&quot;Mercator_1SP&quot;],
PARAMETER[&quot;latitude_of_origin&quot;, 0.0],
PARAMETER[&quot;central_meridian&quot;, 0.0],
PARAMETER[&quot;scale_factor&quot;, 1.0],
PARAMETER[&quot;false_easting&quot;, 0.0],
PARAMETER[&quot;false_northing&quot;, 0.0],
UNIT[&quot;m&quot;, 1.0],
AXIS[&quot;Easting&quot;, EAST],
AXIS[&quot;Northing&quot;, NORTH]]</crs>
<resampling>Nearest</resampling>
<referencePixelX/>
<referencePixelY/>
<easting/>
<northing/>
<orientation/>
<pixelSizeX/>
<pixelSizeY/>
<width/>
<height/>
<tileSizeX/>
<tileSizeY/>
<orthorectify>false</orthorectify>
<elevationModelName/>
<noDataValue>NaN</noDataValue>
<includeTiePointGrids>true</includeTiePointGrids>
<addDeltaBands>false</addDeltaBands>
</parameters>
</node>

<node id="Subset">
<operator>Subset</operator>
<sources>
<sourceProduct refid="Reproject"/>
</sources>
<parameters class="com.bc.ceres.binding.dom.XppDomElement">
<sourceBands/>
<region>0,0,0,0</region>
<referenceBand/>
<geoRegion>POLYGON ((${east} ${north}, ${west} ${north}, ${west} ${south}, ${east} ${south}, ${east} ${north}))</geoRegion>
<subSamplingX>1</subSamplingX>
<subSamplingY>1</subSamplingY>
<fullSwath>false</fullSwath>
<tiePointGridNames/>
<copyMetadata>true</copyMetadata>
</parameters>
</node>

<node id="Write">
<operator>Write</operator>
<sources>
<source>Subset</source>
</sources>
<parameters>
<file>${ofile}</file>
<formatName>NetCDF4-CF</formatName>
</parameters>
</node>

</graph>

------------------------------------------------------

*file.par*

east=-51
west=-63
north=-35
south=-42

referenceComponentPattern=${ORIGINAL_NAME}_chl
dependentComponentPattern=${ORIGINAL_NAME}_sst

Re: mapgen script with MODIS L2

Posted: Mon Jan 30, 2023 5:54 pm America/New_York
by OB SeaDAS - xuanyang02
mapgen on an L2 file is equivalent to l2bin + l3mapgen. It seems mapgen fails doing l2bin with l2_flags as one of the products.

Please try

Code: Select all

python path/SeaDAS/ocssw/bin/mapgen --ifile=path/A2020003181000.L2_LAC_OC.x.nc --ofile=path/mapgenTest_A2020003181000.L2_LAC_OC.x.nc --product=chlor_a --resolution=1km --oformat=netcdf4 --projection=mercator --verbose
The format of of output file is netcdf4.

Re: mapgen script with MODIS L2

Posted: Wed Feb 01, 2023 2:16 pm America/New_York
by fbecker23
Thank you xuanyang02 for your reply.
It is a pity that the flags cannot be processed.

Thank you for your help

Re: mapgen script with MODIS L2

Posted: Thu Feb 02, 2023 9:18 am America/New_York
by OB.DAAC - SeanBailey
The binner cannont bin the flags because they are a bitwise value that cannot be expressed as a floating point value.
The binner uses the l2_flags to exclude data during the binning process - which is what the options you defined in the graph were doing for the band maths step.

The flaguse parameter to l2bin allows the user to set which l2_flags will be used to mask data. The default values are defined to ensure that only the most credible data are used in the L3 output.

The mapgen script is simply a wrapper for running the l2bin and l3mapgen (and if necessary the l2gen) binaries. It was written to produce pretty pictures, so uses a lax set of flags when running l2bin. But you don't need to use the mapgen to do what you seek.

You can accomplish the same exclusion in the graph you defined by passing l2bin:

Code: Select all

flaguse=LAND,CLDICE,CHLFAIL,BOWTIEDEL,STRAYLIGHT

BTW, BOWTIEDEL is only needed for VIIRS, but won't hurt for other sensors.

Then map the resulting bin file with l3mapgen.

If you run mapgen with the -v (verbose) flag you can see how it sets up the calls to l2bin and l3mapgen. You can then copy those and modify the flaguse paramter appropriately in your own workflow.

Re: mapgen script with MODIS L2

Posted: Thu Feb 23, 2023 8:57 am America/New_York
by fbecker23
Excellent!
Thank you Sean!