Hello,

I was just wondering if anyone could clarify what MODIS bands are used in generation of a true color image using msl1tcpcbox by default.  Also, is the result truly a 250m product, or are the other bands (presumably blue and green) interpolated from the 1000 meter bands.  Thank you

Andy

Hi Andy,

It depends on what type of level-1 MODIS file you are starting with.
If you are beginning with a high-res MODIS file, then the default bands
used by msl1tcpcbox are 7, 6, and 2.  Note that msl1tcpcbox plays
fast and loose with band numbers.  For the hi-res MODIS files, see
the following table showing wavelengths and resolutions for each band.

    Pos   Band   Wavelength      Resolution  EV_Position
      0      8       412            1000          0
      1      9       443            1000          1
      2      3       469             500            0
      3     10       488            1000          2
      4     11       531            1000          3
      5     12       551            1000          4
      6      4       555             500            1
      7      1       645             250              0
      8     13       667lo          1000          5
      9     14       678lo          1000          7
     10     15       748            1000          9
     11      2       858             250              1
     12     16       869            1000         10
     13      5      1240             500            2
     14      6      1640             500            3
     15      7      2130             500            4

     16     31      11um            1000
     17     32      12um            1000
     18     22      3.9um           1000
     19     23      4.0um           1000

The "Pos" field gives the msl1tcpcbox band number;
the "Band" field gives the official MODIS band number.

To make a long story short, the default hi-res true color
images use the 250-meter 645-nm band for red, the
500-meter 555-nm band (interpolated to 250-meter) for
green, and the 500-meter 469-nm band (interpolated to
250-meter) for blue.

Regards,

Norman

Hi there,
I notice the msl1tcpcbox program automatically uses the climatology METOZ for atmospheric correction.
How do I get the program to use my METOZ directory which is defined for my login?  (ie if I use MSL12 and define metsearch=1)

Cheers
Matt

Matt,

You can't.  It shouldn't matter much anyway, as it's just a pretty picture.

Sean

You would have to use the l1brsgen command on the Unix command line.  The SeaDAS GUI will not let you change the bands used, but the binary program will allow you to do so.  Assuming you are in the $SEADAS directory, you would execute a command something like this:

./bin/l1brsgen ifile='./demo/A2006302181500.L1B_LAC' geofile='./demo/A2006302181500.GEO' ofile='./A2006302181500.brs' rgb='678,443,412'

Mark

Thanks for the tip!

what's the difference of l1brsgen, l1mapgen vs msl1tcpcbox? I want to generate multiple true color images with projection using commands.

l1brsgen does NOT do projection. 
l1mapgen (formerly msl1tcpcbox) does projection.
rgb rendering is otherwise identical.

The interface for l1mapgen has changed since the previous posts.  See below:

> l1mapgen

l1mapgen 5.9.7 (Apr 22 2009 11:10:28)
Usage: l1mapgen argument-list

The argument-list is a set of keyword=value pairs. The arguments can
be specified on the commandline, or put into a parameter file, or the
two methods can be used together, with commandline over-riding.

The list of valid keywords follows:

   par=input parameter file

   ifile=input L1 file name
   geofile=input L1 geolocation file name (MODIS only)
   ofile=output file name
   resolution=define the processing resolution (MODIS only)
       -1: standard ocean 1km processing, limited band set
     1000: 1km resolution including aggregated 250 and 500m land bands.
      500: 500m resolution including aggregated 250 land bands and replication
          for lower resolution bands.
      250: 250m resolution with replication for lower resolution bands.

   north=northern limit for output image (-90.0 to 90.0)
   south=southern limit for output image (-90.0 to 90.0)
   east=eastern limit for output image (-180.0 to 180.0)
   west=western limit for output image (-180.0 to 180.0)
   width=width of output image in pixels (default=600)
   threshold (0.1 [default])
   rgb=bands to use for red, green and blue (default=sensor specific)
   atmocor=use surface reflectance (default=0, use top of atmosphere reflectance)
   datamin=minumum reflectance for scaling (default=0.01)
   datamax=maximum reflectance for scaling (default=0.9
   stype=set scaling type (default=0)
       0 = log scaled between datamin and datmax
       1 = linearly scaled between datamin and datmax

This program produces a PPM-formatted output image rendered in a
Plate Carree projection.

The default band combination produces a "true-color" image
Other combinations may be chosen with the "rgb=" option.  The expected
argument to this option is a comma separated string of wavelengths that specifies
the desired bands in red-green-blue order.  For example, to produce a false
color SeaWiFS output image using 865, 670 and 555 nm
as the red, green, and blue values respectively, the option would be specified
as "rgb=865,670,555".

Thank you for your quick reply. I have one more question. I read this for generating ture color image.
http://oceancolor.gsfc.nasa.gov/forum/oceancolor/topic_show.pl?tid=2505
The scripts only works for one image scene. I was trying to modify the scripts to do batch processing. I read the scripts. But I don't understand why there is steps 2, 3, and 4 in the scirpts if the command l1mapgen already did the projection. Are these steps only for adding the coast lines and grids and convert to tiff?

Could I omit the steps 2, 3, and 4, and only using the step 1 for the unix command to get true color image with projection?
what projection does "l1mapgen" command used for default? Could I change to a specific projection? Thank you!

One more question is:  I tried to modify the scripts in that link to work for batch processing. but I can't figure out how to use loop to combine both unix commands and seadas commands in the scripts. "seadas -em -b scripts.txt" could not work if I put it in the loop in a unix command file.

Attached are the five steps (from the example script provided in the link http://oceancolor.gsfc.nasa.gov/forum/oceancolor/topic_show.pl?tid=2505):
;(1) Unix command:
cmd1 = 'l1mapgen ifile=A2007244185000.L1B_QKM geofile=A2007244185000.GEO ' +   $
  'ofile=A2007244185000.ppm resolution=250 threshold=0.01 width=3000 ' + $
  'north=50 south=40 west=-95 east=-80 atmocor=1'

; cmd1='msl1tcpcbox -a -- 3000 50 40 -95 -80 0.01 A2007244185000.L1B_QKM A2007244185000.GEO A2007244185000.ppm'
spawn,cmd1
cmd2='pnmtopng A2007244185000.ppm > A2007244185000.png'
spawn,cmd2

;(2) Map L1B file with equatorial cylindrical equidistant
;    projection(default):

load,'A2007244185000.L1B_QKM',gfile='A2007244185000.GEO'
mapimg,band=[1],/CYLIN,LIMIT=[40,-95,50,-80],P0LAT=0.0,xsize=3000, ysize=2000
display,band_no=2,fbuf=1

;(3) Load/Display True-color image. SeaDAS command:
load,'A2007244185000.png',ftype='png'
display,band_no=3,fbuf=2

;(4) Copy navigation info from mapped L1B to True-color image:
;SeaDAS command:

mband_cmd, cmd_array=['result = b3'],navband=2
display,band_no = 4, fbuf=3
coast,color=6
grid, grdcol=5,latdel=2,londel=2

;(5) Output as GeoTiff file:

out,'A2007244185000.tif',ftype='tiff'

One more question. Could the true color image with projection be loaded into Seadas or other software? I hope it can be used to overlay the ground sampling points. Thank you!

Just created a script for your example to handle multiple
image files. The resulted true color png files are portable for
outside SeaDAS use.

Within the For loop, "& $" are needed for statement separator and
line continuation.

Yes, some extra steps were used to copy navigation info for those
mapped images and to create geotiff files. You may skip them if you
can overlay png images on maps.

Long

; *********************************************************************
; This script will create mapped true color images
; from multiple L1B images. And output them as geotiff
; files

files=file_search('A2008*.L1B_QKM',count=n)
print,'No. of files to process:', n

for i=0,n-1 do begin                                                & $
  base_fname = FILE_BASENAME(files(i), '.L1B_QKM')                  & $
  gfile = base_fname + '.GEO'                                       & $
  ppmfile = base_fname + '.ppm'                                     & $
  pngfile = base_fname + '.png'                                     & $
  tiffile = base_fname + '.tif'                                     & $
  ifile   = files(i)                                                & $
  cmd1 = 'l1mapgen ifile=' + ifile + ' geofile=' + gfile +            $
         ' ofile=' + ppmfile  +                                       $
         ' resolution=250 threshold=0.01 width=3000 north=30 ' +      $
         'south=20 west=-85 east=-70 atmocor=1'                     & $
  spawn,[cmd1]                                                      & $
  cmd2='pnmtopng ' + ppmfile + ' > ' + pngfile                      & $
  spawn,[cmd2]                                                      & $
endfor

Thank you for your reply. And thanks for posting the codes. That really helps a lot!
It seems you didn't add the part for generating Tif image.
I tried to add that part before the line endfor, but it doens't work. No tif image were generated. Could you help find what 's wrong?
And also there seems something wrong in these lines:
cmd1 = 'l1mapgen ifile=' + ifile + ' geofile=' + gfile +            $
         ' ofile=' + ppmfile  +                                       $
         ' resolution=250 threshold=0.01 width=3000 north=30 ' +      $

I know there assumed to have + $ for each line. but why there is a + before and after the file name?

;(2) Map L1B file with equatorial cylindrical equidistant projection(default):
load, ifile, gfile=gfile
mapimg,band=[1],/CYLIN,LIMIT=[40,-85,50,-75],P0LAT=0.0,xsize=1000, ysize=1000
display,band_no=2,fbuf=1

;(3) Load/Display True-color image. SeaDAS command:
load, pngfile,ftype='png'
display,band_no=3,fbuf=2

;(4) Copy navigation info from mapped L1B to True-color image:
;SeaDAS command:
mband_cmd, cmd_array=['result = b3'],navband=2
display,band_no = 4, fbuf=3
coast,color=6
grid, grdcol=5,latdel=1,londel=1
;(5) Output as GeoTiff file:

out,tiffile,ftype='tiff'

endfor

A full script including Geotiff output is here:

; This script will create mapped true color images
; from multiple L1B images. And output them as geotiff
; files

files=file_search('A2008*.L1B_QKM',count=n)
print,'No. of files to process:', n

for i=0,n-1 do begin                                                & $
  base_fname = FILE_BASENAME(files(i), '.L1B_QKM')                  & $
  gfile = base_fname + '.GEO'                                       & $
  ppmfile = base_fname + '.ppm'                                     & $
  pngfile = base_fname + '.png'                                     & $
  tiffile = base_fname + '.tif'                                     & $
  ifile   = files(i)                                                & $
  cmd1 = 'l1mapgen ifile=' + ifile + ' geofile=' + gfile +            $
         ' ofile=' + ppmfile  +                                       $
         ' resolution=250 threshold=0.01 width=3000 north=30 ' +      $
         'south=20 west=-85 east=-70 atmocor=1'                     & $
  spawn,[cmd1]                                                      & $
  cmd2='pnmtopng ' + ppmfile + ' > ' + pngfile                      & $
  spawn,[cmd2]                                                      & $
  load,files(i),gfile=gfile                                         & $
  mapimg,band=[1],/CYLIN,LIMIT=[20,-85,30,-70],P0LAT=0.0,xsize=3000, ysize=2000  & $
  display,band_no=2,fbuf=1                                          & $
  load, pngfile,ftype='png'                                         & $
  display,band_no=3,fbuf=2                                          & $
  mband_cmd, cmd_array=['result = b3'],navband=2                    & $
  display,band_no = 4, fbuf=3                                       & $
  coast,color=6                                                     & $
  grid, grdcol=5,latdel=2,londel=2                                  & $
  out,tiffile,ftype='tiff'                                          & $
  clear_up                                                          & $
endfor

;************************************************************************

Due to 24 bit color management limitations in SeaDAS, the 250m high
resolution geotiff image quality is not high. We will improve this soon.

Long   

In SeaDAS command, strings can not mix with variables, so use
" + " to combine them together to a form a command line.
And IDL can "spawn" this command.

Now ifile, gfile and ppmfile are variables in the command line.

cmd1 = 'l1mapgen ifile=' + ifile + ' geofile=' + gfile +            $
         ' ofile=' + ppmfile  +                                     $
         ' resolution=250 threshold=0.01 width=3000 north=30 ' +    $

Long

Thank you. I have generated images sucessfually. But, what I find is that the image in the tif format looks quite different from the png or ppm format. Could you tell me if that is correct? what results in the difference?

Yes, apparently our IDL code does something to the TIFF which reduces it to 256 colors, so it won't look as good as the png or ppm.  We're still working on a fix for this.

Mark

Thank you, Mark. The tiff true color image has the projection. I am wondering if it is possible to project png true color image, so we don't need to use tiff.

No, a PNG is a simple image with no navigation information, so it can't be projected.  Any projections would have to be done prior to creating the PNG.

Mark

Dear Long,

How to execute this full script.
I am not familiar with the scripting.

I named your script as l1mapgeneration. When I called it from linux prompt, I've got the error message as below:

./l1mapgeneration: line 1: syntax error near unexpected token `('
./l1mapgeneration: line 1: `files=file_search('A2010*sub_L1B_QKM',count=n)'

In which folder do I have to put the script in? and from where I have to access? Unix command or seadas command?

Thank you.

Eko

Dear Long,

I now know how to execute the script.
But I had an error message as below.

Using 2130 nm channel for cloud flagging over water.

pnmtopng: Error reading magic number from Netpbm image stream.  Most often, this means your input file is empty.                                                               
Detected MODIS_GDAAC ftype.                                                            
Getting L1B - " EV_250_RefSB_645_01                                                    
" data from HDF file...                                                                
SeaDAS Projection -          128 out of          600 lines processed                   
SeaDAS Projection -          256 out of          600 lines processed                   
SeaDAS Projection -          384 out of          600 lines processed                   
SeaDAS Projection -          512 out of          600 lines processed                   
SeaDAS Projection -          600 out of          600 lines processed                   
% READU: End of file encountered. Unit: 100, File: A2010032033500.png                  
[MAIN]: Error executing the command:                                                   
for i=0,n-1 do begin                                                & base_fname = FILE_BASENAME(files(i), '.L1B_QKM')                  & gfile = base_fname + '.GEO'                                       & ppmfile = base_fname + '.ppm'                                     & pngfile = base_fname + '.png'                                     & tiffile = base_fname + '.tif'                                     & ifile   = files(i)                                                & cmd1 = 'l1mapgen ifile=' + ifile + ' geofile=' + gfile +            ' ofile=' + ppmfile  +                                       ' resolution=500 threshold=0.01 width=600 north=34.0 ' +     ' south=32.5 west=131.0 east=132.5 atmocor=1'              & spawn,[cmd1]                                                      & cmd2='pnmtopng ' + ppmfile + ' > ' + pngfile                      & spawn,[cmd2]                                                      & load,files(i),gfile=gfile                                         & mapimg,band=[1],/CYLIN,LIMIT=[32.5,131.0,34,132.5],P0LAT=0.0,xsize=600, ysize=600  & display,band_no=2,fbuf=1                                          & load, pngfile,ftype='png'                                         & display,band_no=3,fbuf=2                                          & mband_cmd, cmd_array=['result = b3'],navband=2                    & display,band_no = 4, fbuf=3                                       & coast,color=6                                                     & grid, grdcol=5,latdel=2,londel=2                                  & out,tiffile,ftype='tiff'                                          & clear_up                                                          & endfor                                                      
MAIN: Error executing command in batch file

I just changed the long. lat. for my interest region, and 'string' for file_search.
Can you guide me to fix the problem?

Thank you.

Eko

Eko, you can see that pnmtopng encountered an error when creating the png file so please execute each step of the script manually, making sure to check any files created as you go. This will allow you to narrow down the problem.