Re: Extracting a certain data range from a granule of CloudSat data

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Thu May 31 2012 - 14:03:22 MDT

Please, always respond to ncl-talk. THX

==
The most important thing to do in data processing
is to 'look at & understand' the data you will be using.

In NCL, a VERY useful function is "ind"
http://www.ncl.ucar.edu/Document/Functions/Built-in/ind.shtml

====
Data between 20S and 20N only would look like
   http://www.cgd.ucar.edu/~shea/2010128055614_21420_20S-20N.png
Attached:
   2010128055614_21420_CS_2B-GEOPROF_GRANULE_P_R04_E03.tropics.ncl_2
====

If you just want to *directly* read only the data
between your desired latitudes, then you should use the "ind" function
   http://www.ncl.ucar.edu/Document/Functions/Built-in/ind.shtml

Something like

   LAT = eos_file->Latitude_2B_GEOPROF ; (nray_2B_GEOPROF) ==> (37082)
   iTrop = ind(LAT.ge.latS .and. LAT.le.latN) ; desired indicies
   printVarSummary(iTrop) ; size(N)

   data_raw = eos_file->Radar_Reflectivity_2B_GEOPROF(iTrop,:)

and do that fo all appropriate variables. Since the desired subscrips
are not contiguous, this can be quite inefficient. The underlying
HDF software is not particularly efficient at doing this.
(note: the same goes for the underlying netCDF software).

====

My suggestion is to read the entire variable(s) *then* subset
in the main script. This will be much faster.

I can not spend much more time on this.

Good luck

On 5/31/12 9:07 AM, Satomi Sugaya wrote:
>
> Thank you so much your your help.
>
> Yes, I think the image is what I wanted, and the script is really
> helpful. Thank you. I have a further question. In your script, it looks
> like that the tropics range is selected at the stage of plotting. Is
> this true? If it is, I am wondering if you could point me out to a
> function that would let me extract a portion of the whole data for an
> analysis. In other words, from a whole orbit worth of data, I want to be
> able to extract just the tropical part, then save that as a separate
> data set. Can you help me with this?
>
> Thank you so much,
> Satomi
>
> On Wed, May 30, 2012 at 4:58 PM, Dennis Shea <shea@ucar.edu
> <mailto:shea@ucar.edu>> wrote:
>
> Not sure if this is what you want.
>
> See: http://www.cgd.ucar.edu/~shea/__cloudsat_tropics.png
> <http://www.cgd.ucar.edu/~shea/cloudsat_tropics.png>
>
> A sample script is attached.
>
> ===
> Note: that script and the other script uses full 2D vertical
> coordinates.
>
>
> On 5/29/12 4:53 PM, Satomi Sugaya wrote:
>
> Hello,
>
> I am looking at CloudSat data (it's .hdf ), and from one orbit
> worth of
> data, I want to extract one portion (specifically, the portion
> corresponding to tropics) and plot. Following an example
> (provided by
> the HDF group), I have been able to produce a plot with height
> on the
> vertical and the time on the horizontal axes. I want to be able
> to plot
> against latitude and in the tropical regions. I copied my script
> below.
> When I ran it I got the message saying:
>
> ncl 51> plot=gsn_csm_contour(xwks,
> dataf(nbin_2B_GEOPROF|:,nray___2B_GEOPROF|:), res)
> warning:ScalarFieldSetValues: irregular coordinate array sfXArray
> non-monotonic: defaulting sfXArray
> (0)get_lat_values: Warning: Your latitude values do not fall
> between -90 an
> d 90 inclusive.
> (0)You will not get 'nice' latitude labels.
>
> ncl 52>
>
> By testing on the command line, I know that everything worked
> until the
> very last line, where I try to make the plot.
>
> Also, I can't find a reference page for specifying res@...
> variables.
> Could you please point me out to this also?
>
> Thank you so much,
> Satomi
>
> PS The data used here can be downloaded from:
> http://www.hdfeos.org/zoo/__index_openCDPC_Examples.php#__table2
> <http://www.hdfeos.org/zoo/index_openCDPC_Examples.php#table2>
> it's the CloudSat swath data.
>
>
> What I did:
>
> load "$NCARG_ROOT/lib/ncarg/nclex/__gsun/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/__nclscripts/csm/gsn_csm.ncl"
> load "$NCARG_ROOT/lib/ncarg/__nclscripts/csm/contributed.__ncl"
>
>
> begin
> file_name =
> "2010128055614_21420_CS_2B-__GEOPROF_GRANULE_P_R04_E03.hdf"
>
> eos_file=addfile(file_name+ ".he2","r")
>
> data_raw = eos_file->Radar_Reflectivity___2B_GEOPROF
>
> lat = eos_file->Latitude_2B_GEOPROF
> lat@long_name = "Latitude"
> lat@units = eos_file@Latitude_units_2B___GEOPROF
>
> lon = eos_file->Longitude_2B_GEOPROF
>
> lev = eos_file->Height_2B_GEOPROF
>
> lev@long_name = "Height"
> lev@units = eos_file@Height_units_2B___GEOPROF
>
> time = eos_file->Profile_time_2B___GEOPROF
> time@long_name = eos_file@Profile_time_long___name_2B_GEOPROF
> time@units = eos_file@Profile_time_units___2B_GEOPROF
>
> hdf4_file=addfile(file_name,"__r")
>
> data_hdf4 = hdf4_file->Radar_Reflectivity
>
> dataf = tofloat(data_raw)
> dataf = dataf / data_hdf4@_FillValue
> dataf@_FillValue = data_hdf4@_FillValue
>
> dataf = where(dataf.gt.data_hdf4@__valid_range(0) .and.
> dataf.lt.2000,
> dataf, data_hdf4@_FillValue)
>
> dataf!0 = data_raw!0
> dataf!1 = data_raw!1
>
> dataf&nbin_2B_GEOPROF = lev(0,:)
> dataf&nray_2B_GEOPROF = lat
>
> xwks = gsn_open_wks("pdf", "test.ncl") ; open workstation
>
> gsn_define_colormap(xwks,"__BlAqGrYeOrReVi200") ; define colormap
>
> res = True ; plot mods desired
> res@cnFillOn = True ; enable contour fill
> res@gsnMaximize = True ; make plot large
> res@gsnPaperOrientation = "portrait" ; force portrait orientation
> res@cnLinesOn = False ; turn off contour lines
> res@cnLineLabelsOn = False ; turn off contour line
> labels
> res@gsnSpreadColors = True ; use the entire color
> spectrum
> res@cnFillMode = "RasterFill" ; faster
> res@lbOrientation = "vertical" ; vertical labels
> res@cnMissingValFillPattern = 0 ; missing value pattern is set to
> "SolidFill"
> res@cnMissingValFillColor = 0 ; white color for missing values
> res@lbLabelAutoStride = True ; ensure no label overlap
> res@tiMainString = file_name
> res@tiXAxisString = lat@long_name+" ("+lat@units+")"
> res@tiYAxisString = lev@long_name+" ("+lev@units+")"
> res@gsnLeftString=data_hdf4@__long_name
> res@gsnRightString=data_hdf4@__units
>
> plot=gsn_csm_contour(xwks,
> dataf(nbin_2B_GEOPROF|:,nray___2B_GEOPROF|:), res)
>
> end
>
>
> _________________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/__mailman/listinfo/ncl-talk
> <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
>
>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk

Received on Thu May 31 14:03:39 2012

This archive was generated by hypermail 2.1.8 : Wed Jun 06 2012 - 15:17:44 MDT