plotting CALIPSO data

Hi
I am getting into all sorts of trouble plotting the CALIPSO VFM data. The
arrays are very large 4400*5000, and I get error messages:
(0) 0 CAL_LID_L2_VFM-Prov-V2-01.2006-06-13T07-20-12ZD.hdf
warning:ScalarFieldSetValues: coordinate array sfXArray requires 24964297
elements: defaulting
warning:ScalarFieldSetValues: coordinate array sfYArray requires 24964218
elements: defaulting

How do we just subset the plot into small region based on lat and long, so
that I can plot it easily.
The lat , long are in separate arrays and are 1-D, so I could not use
closest_val to get the location of indices from where to subset as well.

Looking forward for you kind suggestion.
Regards,
Prabhakar

PS: I have also attached the code I am using.

;hdfreadCALIOP.ncl
;HERE WE READ AND PLOT THE VFM V2.0 products

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"

begin

 diri = "/Data/CALIPSO/"
 fili = systemfunc("cd "+diri+" ;ls CAL* " )
 nfil = dimsizes( fili )

 do nf = 0,nfil-1 ;BIG LOOP
 print(nf+" "+fili(nf))
 f = addfile(diri+fili(nf), "r")

 var = f->Feature_Classification_Flags

 N = dimsizes(var)
 height = ispan(0,1000,1)

; GEOLOCATION ARRAYS

 lat2d = f->Latitude(fakeDim0|:,fakeDim1|0)
 lon2d = f->Longitude(fakeDim2|:,fakeDim3|0)

; PLOTTING

  if (nf.lt.10) then
   wks = gsn_open_wks ("eps", "vfm_00"+sprintf("%1.0f",nf))
  else
   if (nf.lt.100) then
    wks = gsn_open_wks ("eps", "vfm_0"+sprintf("%2.0f",nf))
   else
    wks = gsn_open_wks ("eps", "vfm_"+sprintf("%3.0f",nf))
   end if
  end if

; colors =
(/"red","RoyalBlue","LightSkyBlue","darkorange","yellow","yellowgreen","thistle","black"/)
; gsn_define_colormap(wks,colors)

; START PANEL

 plot = new(2,graphic)

 res = True

 res@gsnDraw = False ; don't draw
 res@gsnFrame = False ; don't advance frame
 res@gsnMaximize = True ; maximize pot in frame

; SPECIFY SUBPLOT SPECIFIC GEOLOCATION AND TITLE
 res1 = res
 res1@tiMainString = "VFM"

 res1@mpOutlineBoundarySets = "National"
 res1@pmTickMarkDisplayMode = "Always"
 res1@mpDataBaseVersion = "MediumRes" ; Higher res coastline
 res1@mpLimitMode = "LatLon"
 res1@mpMinLatF = 20.
 res1@mpMaxLatF = 40.
 res1@mpMinLonF = 60.
 res1@mpMaxLonF = 90.

;
 res2 = res
 res2@tiMainString = "VFM"
 res2@cnFillOn = True ; color Fill
 res2@cnFillMode = "RasterFill" ; Raster Mode
 res2@cnRasterSmoothingOn = True
; res2@cnRasterMinCellSizeF = 0.0005
 res2@cnLinesOn = False ; Turn off contour lines
 res2@cnLineLabelsOn = False ; Turn off contour lines

 res2@sfXArray = lat2d
 res2@sfYArray = height
 res2@trGridType = "TriangularMesh"
 res2@gsnAddCyclic = False ; Data is not cyclic
;
 res2@lbLabelStrings = ispan(0,7,1) ; 0, 1, 2, 3, 4
 res2@lbLabelPosition = "Center" ; label position
 res2@lbLabelAlignment = "BoxCenters"
 res2@lbTitleString =" "; "0=low/no confidence, 1=clear air, 2=cloud ,
3=aerosol, 4=stratospheric features,/
                             ; 5=surface, 6=subsurface, 7=no signal"
 res2@lbTitlePosition = "Bottom"
 res2@lbTitleFontHeightF = 0.0125 ; default 0.025

;
 plot(0) = gsn_csm_map_ce(wks,res1) ; create map
 pres = True
 pres@gsLineThicknessF = 3.0
 pres@gsLineColor ="red"
 line = gsn_add_polyline(wks,plot(0),lon2d,lat2d,pres)

 plot(1) = gsn_csm_contour(wks,var, res2)

 resP = True
 resP@txString = date_str
 gsn_panel(wks,plot,(/2,1/),resP)
 draw(plot)
 frame(wks)
 end do
end

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Sun May 2 13:01:04 2010