create_fake_mod06.pro

IDL源码

PRO CREATE_FAKE_MOD06, MOD021KM, MOD06, OUTFILE

;+
; Creates a fake MOD06 HDF file that can fool the McIDAS MOD35 ADDE server.
; The ADDE server reads only the following SDS arrays from a MOD06 HDF file:
; 'Latitude', 'Longitude', 'Surface_Temperature'
;
; MOD021KM    Name of input MOD021KM HDF file (from DAAC or IMAPP)
; MOD06       Name of input MOD35 binary file (from DB version of cloud mask)
; OUTFILE     Name of output MOD35 HDF file (fake version)
;-

;- Check arguments
if (n_elements(mod021km) eq 0) then message, 'Argument MOD021KM is undefined'
if (n_elements(mod06) eq 0)    then message, 'Argument MOD06 is undefined' 
if (n_elements(outfile) eq 0)  then message, 'Argument OUTFILE is undefined' 
;-------------------------------------------------------------------------------
; READ INPUT
;-------------------------------------------------------------------------------

;- Read the lat/lon data from the MOD021KM file
hdfid = hdf_sd_start(mod021km)
hdf_sd_varread, hdfid, 'Latitude', lat
hdf_sd_varread, hdfid, 'Longitude', lon
hdf_sd_end, hdfid

;- Remove the last value from each row
lat = lat[0:269, *]
lon = lon[0:269, *]


;- Get the number of pixels and lines
dims = size(lat, /dimensions)
nx = 1354
ny = dims[1] * 5L

;- Get the cloud top properties data from the flat file
openr, lun, mod06, /get_lun
ctp = fltarr(nx/5, 48, ny/5)
readu, lun, ctp
free_lun, lun

bt=intarr(nx/5,7,ny/5)
for i = 0 , 6, 1 do begin
   bt(*,i,*) = round(ctp(*,i,*) * 100. - 15000. )
endfor
;- Rearrange Brightness Temperature array 
brite=transpose(bt,[0,2,1])

sfctemp=intarr(nx/5,1,ny/5)
sfctemp=round(ctp(*,7,*) * 100. - 15000. )
sfct=fix(reform(sfctemp,nx/5,ny/5))

sfcpres=intarr(nx/5,1,ny/5)
sfcpres=round(ctp(*,8,*) * 10.)
sfcp=fix(reform(sfcpres,nx/5,ny/5))
loc=where(sfcp lt 0, count)
if count gt 0 then sfcp(loc) = -32768

pf=intarr(nx/5,1,ny/5)
pf=fix(ctp(*,9,*))
pflag=byte(reform(pf,nx/5,ny/5))
loc=where(pflag lt 0, count)
if count gt 0 then pflag(loc) = 127

cm=intarr(nx/5,1,ny/5)
cm=fix(ctp(*,10,*))
cmethod=byte(reform(cm,nx/5,ny/5))
loc=where(cmethod gt 6, count)
if count gt 0 then cmethod(loc) = 127

ctp1=intarr(nx/5,1,ny/5)
ctp1=round( ctp(*,11,*) * 10. )
cloudtop=fix(reform(ctp1,nx/5,ny/5))
loc=where(cloudtop lt 0, count)
if count gt 0 then cloudtop(loc) = -32768

ctp3=intarr(nx/5,1,ny/5)
ctp3=round( ctp(*,12,*) * 10. )
cloudtopnight=fix(reform(ctp3,nx/5,ny/5))
loc=where(cloudtopnight lt 0, count)
if count gt 0 then cloudtopnight(loc) = -32768

ctp2=intarr(nx/5,1,ny/5)
ctp2=round( ctp(*,13,*) * 10. )
cloudtopday=fix(reform(ctp2,nx/5,ny/5))
loc=where(cloudtopday lt 0, count)
if count gt 0 then cloudtopday(loc) = -32768

ctt1=intarr(nx/5,1,ny/5)
ctt1=round( ctp(*,14,*) * 100. - 15000. )
ctemp=fix(reform(ctt1,nx/5,ny/5))
loc=where(ctemp lt 0 or ctemp eq 17768, count)
if count gt 0 then ctemp(loc) = -32768

ctt2=intarr(nx/5,1,ny/5)
ctt2=round( ctp(*,15,*) * 100. - 15000. )
ctempnight=fix(reform(ctt2,nx/5,ny/5))
loc=where(ctempnight lt 0 or ctempnight eq 17768, count)
if count gt 0 then ctempnight(loc) = -32768

ctt3=intarr(nx/5,1,ny/5)
ctt3=round( ctp(*,16,*) * 100. - 15000. )
ctempday=fix(reform(ctt3,nx/5,ny/5))
loc=where(ctempday lt 0 or ctempday eq 17768, count)
if count gt 0 then ctempday(loc) = -32768

th=intarr(nx/5,1,ny/5)
th=round(ctp(*,17,*) * 10.)
trop=fix(reform(th,nx/5,ny/5))
loc=where(trop lt 0, count)
if count gt 0 then trop(loc) = -32768

cf1=intarr(nx/5,1,ny/5)
cf1=fix(ctp(*,18,*) * 100. )
cfrac1=byte(reform(cf1,nx/5,ny/5))
loc=where(cfrac1 gt 100, count)
if count gt 0 then cfrac1(loc) = 127

cf2=intarr(nx/5,1,ny/5)
cf2=fix(ctp(*,19,*) * 100. )
cfracnight=byte(reform(cf2,nx/5,ny/5))
loc=where(cfracnight gt 100, count)
if count gt 0 then cfracnight(loc) = 127

cf3=intarr(nx/5,1,ny/5)
cf3=fix(ctp(*,20,*) * 100. )
cfracday=byte(reform(cf3,nx/5,ny/5))
loc=where(cfracday gt 100, count)
if count gt 0 then cfracday(loc) = 127

ce1=intarr(nx/5,1,ny/5)
ce1=fix(ctp(*,21,*) * 100. )
cee1=byte(reform(ce1,nx/5,ny/5))
loc=where(cee1 lt 1, count)
if count gt 0 then cee1(loc) = 127

ce2=intarr(nx/5,1,ny/5)
ce2=fix(ctp(*,22,*) * 100. )
cee2=byte(reform(ce2,nx/5,ny/5))
loc=where(cee2 lt 1, count)
if count gt 0 then cee2(loc) = 127

ce3=intarr(nx/5,1,ny/5)
ce3=fix(ctp(*,23,*) * 100. )
cee3=byte(reform(ce3,nx/5,ny/5))
loc=where(cee3 lt 1, count)
if count gt 0 then cee3(loc) = 127

ctp4=intarr(nx/5,1,ny/5)
ctp4=round( ctp(*,24,*) * 10. )
ctpir=fix(reform(ctp4,nx/5,ny/5))
loc=where(ctpir lt 0, count)
if count gt 0 then ctpir(loc) = -32768

scf=intarr(nx/5,5,ny/5)
for i = 0 , 4, 1 do begin
   j = 25 + i 
   scf(*,i,*) = round(ctp(*,j,*) * 100. )
endfor
;- Rearrange Spectral Cloud Forcing Array
spectral=transpose(scf,[0,2,1])

ctp5=intarr(nx/5,5,ny/5)
for i = 0 , 4, 1 do begin
   j=30+i
   ctp5(*,i,*) = round(ctp(*,j,*) * 10. )
endfor
;- Rearrange Cloud_Top_Pressure_From_Ratios Array
ctpratios=transpose(ctp5,[0,2,1])

sfc=intarr(nx/5,1,ny/5)
sfc=round( ctp(*,35,*) )
sfctype=fix(reform(sfc,nx/5,ny/5))
loc=where(sfctype lt 0, count)
if count gt 0 then sfctype(loc) = -32768

rv=intarr(nx/5,7,ny/5)
for i = 0, 6, 1 do begin
   j = 36 + i 
   rv(*,i,*) = round(ctp(*,j,*) * 100. )
endfor
;- Rearrange Radiance Variance Array
variance=transpose(rv,[0,2,1])

btdif=intarr(nx/5,2,ny/5)
for i = 0, 1, 1 do begin
   j = 43 + i
   btdif(*,i,*) = round(ctp(*,j,*) * 100. )
endfor
;- Rearrange Brightness Temperature Difference Array
difference=transpose(btdif,[0,2,1])

ph1=intarr(nx/5,1,ny/5)
ph1=fix(ctp(*,45,*) )
phase1=byte(reform(ph1,nx/5,ny/5))
loc=where(phase1 gt 100, count)
if count gt 0 then phase1(loc) = 127

ph2=intarr(nx/5,1,ny/5)
ph2=fix(ctp(*,46,*) )
phase2=byte(reform(ph2,nx/5,ny/5))
loc=where(phase2 gt 100, count)
if count gt 0 then phase2(loc) = 127

ph3=intarr(nx/5,1,ny/5)
ph3=fix(ctp(*,47,*) )
phase3=byte(reform(ph3,nx/5,ny/5))
loc=where(phase3 gt 100, count)
if count gt 0 then phase3(loc) = 127

;-------------------------------------------------------------------------------
; WRITE OUTPUT
;-------------------------------------------------------------------------------

;- Open the output HDF file
hdfid = hdf_sd_start(outfile, /create)

;- Write the latitude
varid = hdf_sd_create(hdfid, 'Latitude', size(lat, /dimensions), /float)
hdf_sd_adddata, varid, lat
hdf_sd_endaccess, varid

;- Write the longitude
varid = hdf_sd_create(hdfid, 'Longitude', size(lon, /dimensions), /float)
hdf_sd_adddata, varid, lon
hdf_sd_endaccess, varid


;- Write the cloud top properties SDS's  
varid = hdf_sd_create(hdfid, 'Brightness_Temperature', size(brite, /dimensions), /int)
hdf_sd_adddata, varid, brite
hdf_sd_attrset, varid, 'units', 'K'
hdf_sd_attrset, varid, 'scale_factor', .01d0, /double 
hdf_sd_attrset, varid, 'add_offset', -15000., /double 
hdf_sd_attrset, varid, 'valid_range', [0, 20000], /short
hdf_sd_attrset, varid, '_FillValue', -32768, /short
hdf_sd_endaccess, varid