Not logged inOcean Color Forum
Forum Ocean Color Home Help Search Login
Up Topic SeaDAS / SeaDAS 6.x - General Questions / smigen (locked)
- By c.binding Date 2007-11-14 16:28
I have just tested the smigen function below, which works fine but for the fact that the no-data pixels that are zeros in the original L3b file are re-assigned the value of datamax in the smi.  How can I force them to remain zeros?

smigen ifile='A2005145.daily_L3b_2km' prod='nLw_551' meas=1 ofile='A2005145.daily_L3b_2km.smi' stype=1 projection=RECT resolution=2km precision=I latnorth=50.00 latsouth=40.00 lonwest=-94.00 loneast=-74.00 datamin=0 datamax=5 gap_fill=2

Caren
- By long Date 2007-11-14 17:17
Caren,

This is designed to re-assign no-data values from 0 to 65535 in SMIGEN.
Users has no controls for those data values at program run.

However, users always have freedom to reset those value as 0 (or any other values)
with a few SeaDAS commands:

Main Menu-> Utilities -> Data Manipulation -> User Defined Operation

In command window, just type in commands to create a new image:

  result = B1                                            ; copy band 1
  ix = WHERE( B1 EQ 65535, icnt)            ; find data value=65535           
  if (icnt GT 0) then result(ix) = 0.0             ; reset those 65535 values as 0.0

You may also refer other examples provided within this function for more
advanced calculations.

Long
- By mike Date 2007-11-14 17:25
And for automation you could use the mband_cmd command (see the FAQ for an example).
- By c.binding Date 2007-11-15 09:33
OK, thanks Long.
- By c.binding Date 2007-11-15 14:00
I've incorporated smigen into a batch program that ends with plotting, annotating and export the smi to .png.

the smi is produced correctly but the .png is produced with a global coastline overlaid.  If I do the annotations on the same smi from the GUI then it works fine.  the section of script is below - am I missing something with loading/displaying the smi?

spawn, 'smigen ifile=' + l3b_file + ' prod=nLw_551 meas=1 ofile=' + smi_file + ' stype=1 projection=RECT resolution=2km precision=I latnorth=50.00 latsouth=40.00 lonwest=-94.00 loneast=-74.00 datamin=0 datamax=100 gap_fill=2'

seadisp
load, smi_file
mband_cmd, navband=1, cmd_array=['result=B1', 'ix=WHERE(B1 GT 100, icnt)', 'if (icnt GT 0) then result(ix) = 0.0']

  loadpal, '/media/IOMEGA_HDD1/modis/nlw_scale'
  display, band_no=2, fbuf=1, smin=0, smax=5
  landmask, color=7
  coast, color=6
  grid, grdcol=6, glinestyle=2, label=0
  title = date
  text ='MODIS L3b nLw551'
  outfile = smi_file + '.png'
  xyouts,10,10,title,color=16,alignment=0,charsize=2,charthick=2,/device
  xyouts,10,50,text,color=16,alignment=0,charsize=2,charthick=2,/device
  tvlct, r, g, b, /get
  img_new = TVRD(true=1)
  write_png, outfile, img_new, r, g, b
  clear_up
- By long Date 2007-11-15 16:57 Edited 2007-11-15 17:13
Your script seems fine to me, because I used it with SeaDAS 5.1
and able to produced a good png picture.

Can you do a test with my script (with a little change) to see any
better results ? Maybe a newer SeaDAS is the solution.
SeaDAS 5.1 does not use 'seadisp' any more.

Below is the working script:

load,'A20072652007272.L3m_8D_CHLO_4',ftype='modis', limit=[40,-94,50,-74], prod_ix=[1]
mband_cmd, navband=1, cmd_array=['result=B1', 'ix=WHERE(B1 GT 100, icnt)', 'if (icnt GT 0) then result(ix) = 0.0']
display
loadpal,'/home/lwang/seadas5.1/data/common/luts/standard/02-standard_chl.lut'
landmask, color=7
coast, color=6
grid, grdcol=6, glinestyle=2, label=0
title = 'Oct. 15, 2007'
text ='MODIS L3b nLw551'
outfile = 'test.png'
xyouts,10,10,title,color=16,alignment=0,charsize=2,charthick=2,/device
xyouts,10,50,text,color=16,alignment=0,charsize=2,charthick=2,/device
tvlct, r, g, b, /get
img_new = TVRD(true=1)
write_png, outfile, img_new, r, g, b
clear_up
Attachment: test.png (16k)
- By c.binding Date 2007-11-16 09:43
ok, I tried your script and got even more bizzare results.

I've attached images - test from the script you sent and test_GUI using the same smi.

I'll try updating SeaDAS - does the new version still have the gap_fill function in smigen?
Attachment: test_GUI.png (32k)
Attachment: test.png (550B)
- By mike Date 2007-11-16 17:05
Yes the gap_fill option will still be available to you. Also, if you happen to have a full IDL license, on the SeaDAS command-line try manually executing each of the script's commands one-by-one in order to see if you get the same problem.
- By c.binding Date 2007-11-19 15:01
script works fine on seadas 5.1

thanks
- By long Date 2007-11-19 18:18
Good to know that SeaDAS 5.1 works for you better.
- By c.binding Date 2007-11-20 15:46
I lied.

I have a really strange problem with the script now.  my script takes l3b images, creates a daily binned file, applies smigen and then plots the smi with annotations and exports as png.  The script stalls loading the smi with the error message:

Output File                     : /media/SDB1/MODIS/2005/May/A2005145.daily_smi
Detected MODIS ftype.
% HDF_SD_START: Invalid HDF file or filename ((null)).

The problem appears to be in defining the paths for l3b_file and smi_file (bold text in the script below), particularly with the parameter MMM.

if I replace the MMM with, for example, the string 'May' then the script works fine and creates the png file.  Using MMM, l3b_file and smi_file are created, but it stalls loading smi_file.

I'm sure it'll be something simple but I've exhausted all my troubleshooting ideas.  Your help would be much appreciated again!

months = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']
dir='/media/IOMEGA_HDD/modis/downloads/'
dir_out='/media/SDB1/MODIS/'

bs=strlen(dir)

filenames = findfile(dir + '*.L3b_2km')
n=n_elements(filenames)

start_day=uint(strmid(filenames(0),bs+5,3))
end_day=uint(strmid(filenames(n-1),bs+5,3))

for j=start_day, end_day do begin & $
jday=strtrim(j,1) & $
dl=uint(strlen(jday)) & $
if (dl EQ 1) then begin & $
  jday=('00'+jday) & $
endif else if (dl EQ 2) then begin & $
  jday=('0'+jday) & $
endif else begin & $
  jday=jday & $
endelse & $

   ;***Extracting L3 files for specified day and creating filelist***
for i=0,n-1 do begin & $
   DDD = strmid(filenames(i),bs+5,3) & $
   if (DDD EQ jday) THEN BEGIN & $
     cmd = 'echo ' + filenames(i) + ' >>filelist.txt' & $
     SPAWN,[cmd] & $
     ;***Creating DAY & MONTH strings***
     YYYY = strmid(filenames(i),bs+1,4) & $
     DDD = strmid(filenames(i),bs+5,3) & $
     cmd = '$SEADAS/bin/jd2date' + ' ' + DDD + ' ' + YYYY & $
     spawn, cmd, datestring & $
     DAY = strmid(datestring,2,2) & $
     MONTH = strmid(datestring,0,2) & $
     MMM = months(fix(MONTH)-1) & $
     date = DAY + ' ' + MMM + ' ' + YYYY & $
   endif & $
endfor & $

spawn,'ls -l '+'filelist.txt',output & $
IF output[0] EQ '' THEN BEGIN & $
   print, 'There is no imagery for   ', title, '.  Skipping to next day.' & $
   goto, ioerr & $
ENDIF & $

l3b_file=(dir_out+YYYY+'/' + MMM + '/A' + YYYY + jday + '.daily_L3b_2km') & $
smi_file=(dir_out+YYYY+'/' + MMM + '/A' + YYYY + jday + '.daily_smi') & $


l3bin, ifile='filelist.txt', ofile=l3b_file, l3bprod=  ['nLw_412', 'nLw_443', 'nLw_488', 'nLw_531', 'nLw_551', 'nLw_667', 'nLw_678', 'nLw_748', 'nLw_869'], noext=1 & $

spawn, 'smigen ifile=' + l3b_file + ' prod=nLw_551 meas=1 ofile=' + smi_file + ' stype=1 projection=RECT resolution=2km precision=I latnorth=50.00 latsouth=40.00 lonwest=-94.00 loneast=-74.00 datamin=0 datamax=100 gap_fill=2' & $

load,smi_file & $
mband_cmd, navband=1, cmd_array=['result=B1', 'ix=WHERE(B1 GT 100, icnt)', 'if (icnt GT 0) then result(ix) = 0.0'] & $
display & $
loadpal,'/media/IOMEGA_HDD/modis/color_scales/nlw_scale' & $
landmask, color=7 & $
coast, color=6 & $
grid, grdcol=6, glinestyle=2, label=0 & $
title = date & $
title2='Environment Canada' & $
title3='Remote Sensing Group' & $
text ='MODIS Daily nLw551' & $
outfile = 'test.png' & $
xyouts,20,40,title,color=16,alignment=0,charsize=1.5,charthick=1.5,/device & $
xyouts,750,40,title2,color=16,alignment=0,charsize=1.5,charthick=1.5,/device & $
xyouts,750,10,title3,color=16,alignment=0,charsize=1.5,charthick=1.5,/device & $
xyouts,20,10,text,color=16,alignment=0,charsize=1.5,charthick=1.5,/device & $
tvlct, r, g, b, /get & $
img_new = TVRD(true=1) & $
write_png, outfile, img_new, r, g, b & $
clear_up & $

SPAWN,['rm -f filelist.txt'] & $

ioerr:  & $

endfor
- By long Date 2007-11-20 17:56
Caren,

I am impressed by your nicely written command script. This is really
good for data processing automation. I learn from your nice code.

One suggestion may be helpful to self troubleshoot/debug is use of
print statement between command lines to monitor the processing
steps (dump variable values with screen messages).

I created some dummy l3b files and run your first part of script with
some "print statement" help, those input/output files and paths can
be checked step by step.

Long

months = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']
dir='/home/long/test/'
dir_out='/home/long/modis/'

bs=strlen(dir)

filenames = findfile(dir + '*.l3b')

n=n_elements(filenames)

print,'N files = ', n
print,'files = ', filenames

start_day=uint(strmid(filenames(0),bs+5,3))
end_day=uint(strmid(filenames(n-1),bs+5,3))

for j=start_day, end_day do begin & $
  jday=strtrim(j,1) & $
  . . . . .

  l3b_file=(dir_out+YYYY+'/' + MMM + '/A' + YYYY + jday + '.daily_L3b_2km') & $
  smi_file=(dir_out+YYYY+'/' + MMM + '/A' + YYYY + jday + '.daily_smi') & $

  print,'smi_file = ', smi_file & $
  . . . . . .

Screen Messages:N files =            3
files =  /home/long/test/A20073072007309.l3b /home/long/test/A20073072007310.l3b
/home/long/test/A20073072007311.l3b
smi_file =  /home/long/modis/2007/Nov/A2007307.daily_smi
. . . . .

- By haidichen Date 2009-11-14 10:32
Hi, everyone, I am vary appreciate that someone could help me to see this question as it is a bit confusion.
1. first, for L3bin, I binned chlor_a data into three days composite, 1km. Suppose I binned the day 1,2,3. If only day 1 and 3 have meaningful value in that pixel,  so the result will be [chlor_a(1) +chlor_a(2) +chlor_a(3)]/weight, generally the weight is sqrt(N), N=days, so here N=3, then won't the missing value in chlor_a(2) influence the result?

2. Later I use smigen to output the data into SMI format, and I try two different projection: SIN and RECT. the result is different and SIN has more pixels and can not match with the coast line and land mask. So, could anyone tell me what is the difference between SIN and RECT? Which has a higher quality?

Thank you so much if any one can help me.
- By mike Date 2009-11-16 11:20
1. N does not equal the number of days. It equals the number of valid pixels for that bin. So for a missing chlor_a(2) bin, the result will be [chlor_a(1) + chlor_a(3)]/sqrt(2).

2. Can you please email a small sample SMI file where this problem occurs and tell us the exact SIN and RECT projection parameters you are using?:  seadas [at] seadas.gsfc.nasa.gov
- By Caiyun Date 2010-03-11 15:26
I have problem using smigen. My par file is:

ifile=2002335072921.N16.l2bin.avhrrsst.hdf
ofile=2002335072921.N16.smigen.avhrrsst.hdf
prod=avhrr_sst
precision=I
datamin=-2.0
datamax=45.0
stype=1
meas=1
loneast=180.000
lonwest=-180.000
latnorth=90.0000
latsouth=-90.0000
projection=RECT
resolution=4km
proddesc="avhrr_sst"
units="Deg C"

I always get "cannot find product name", and cannot locate palette

Please help me.

Caiyun
- By long Date 2010-03-11 16:51
Can you provide your 2002335072921.N16.l2bin.avhrrsst.hdf file
for a test please ? You may compress it and send in email or
in attachment.

Long
Up Topic SeaDAS / SeaDAS 6.x - General Questions / smigen (locked)



Responsible NASA Official: Gene C. Feldman
Curator: OceanColor Webmaster
Authorized by: Gene C. Feldman
Updated: 03 July 2013
Privacy Policy and Important Notices NASA logo