Fwd: Re: passing arguments, contouring

From: Jeremy Winick <winickjr_at_nyahnyahspammersnyahnyah>
Date: Fri Jun 14 2013 - 16:48:08 MDT

Thank you for the information. There are many ways of "skinning the cat"
and as I work more with NCL and get out of the Matlab think, I will feel
more comfortable. The use of str_get_field seems a way to get the one
field (column) of a 24 column long (60000 line) data file that isn't
numeric, but then using it many more times to get other fields then
string to float seems inefficient. The readAsciiTable gets everything
except the one string field (which happens to be a month-day
abbreviation (such as Feb12)). Could I use the str_get_field to extract
the one column, then use readAsciiTable for the rest? Does a file need
to be closed or rewound in an NCL script to do this?

I will try the contour plot with the one-dim temperature field and the
sfYArray and sFXarray to the lat and lon variables. Do these still need
the units "degrees_north" and "degrees_east" ? The first attempt this
way didn't work, I get:
(0) Error: scalar_field: If the input data is 1-dimensional, you must
set sfXArray and sfYArray to 1-dimensional arrays of the same length.
warning:create: Bad HLU id passed to create, ignoring it

As far as I can tell, I have set the sfYArray and sfXArray, but I don't
understand the warning about bad HLU id passed. The script seems to be
hanging after this error message, and I hit control-C to kill it.

Here is the code, with portions that set up the grid from natgrid
commented out.

; first test ncl script to read in ascii SABER data, say
; then create ncl variables for temperatures at various pressure levels,
associated lat and lon
; suggestion in this script modify so not using natgrip, but use
sfYArray set to original lat points and sfXArray set to orig lon pts.
; then try contour plot on polar stereographic map
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"
; the contour map routines such as gsn_csm_contour_map_polar seem to
want data to be 2-D and have 1-D lat and long coordinate
; arrays attached, and longitude coord array must have valid "units"
attribute of something like "degrees_east" (similarly for lat)
; lat array same length as leftmost dimension of the 2-D array and
longitude array must be same as the rightmost dim of data (2-d array).
; thus raw satellite data does not do this- it needs to be gridded
; 9 June 2013 works for hand-wired separating out first day. Next choose
day or days by input parameter. Choose level by input parameter?
: There is an ncl function day_of_year(year:integer, month:integer,
   system("cd diri")
   data =
   latd = data(:,5)
; print(lat)
   lond = data(:,6)
   doy = data(:,2)
   t50mb = data(:,9)
   t20mb = data(:,10)
   t10mb = data(:,11)
   t5mb = data(:,12)
; test for essentially the first day, 1-477 profiles
   t1 = t50mb(0:476)
   lat1 = latd(0:476)
   lon1 = lond(0:476)
   do n=0,476
    if(lon1(n).lt.0.0) then
       lon1(n) = 360.0 + lon1(n)
   end if
   end do
; create lat-lon grid
; latp= fspan(30, 85, 12)
; printVarSummary(latp)
; lonp= fspan(0,360,37)
; use natgrid to grid original satellite data
; tgrid=natgrid(lat1,lon1,t1,latp, lonp)
; write_matrix(tgrid,"12f7.2", False)
; print(latp)
; print(lonp)
; printVarSummary(tgrid)
; printVarSummary(latp)
; printVarSummary(lonp)
; tgrid!0 = "lat"
; tgrid!1 = "lon"
; tgrid@units='Degrees K"
   lat1@units = "degrees_north"
   lon1@units = "degrees_east"

; tgrid&lat = latp
; tgrid&lon = lonp
   sfXArray = lon1
   sfYArray = lat1
  wks = gsn_open_wks("x11", "polar")
  res = True
  res@gsnAddCyclic = False
  res@gsnPolar = "NH"
  res@mpMinLatF = 30
  res@cnFillOn = True

plot = gsn_csm_contour_map_polar(wks,t1,res)
   wks = gsn_open_wks("ps", "polar")
  res = True
  res@gsnAddCyclic = False
  res@gsnPolar = "NH"
  res@mpMinLatF = 30
  res@cnFillOn = True

plot = gsn_csm_contour_map_polar(wks,t1,res)
   wks = gsn_open_wks("png", "polar")
; note in above wks open "polar" is file name for output , this case
would be polar.png
; make name specific to case so each new plot has distinct name
  res = True
  res@gsnAddCyclic = False
  res@gsnPolar = "NH"
  res@mpMinLatF = 30
  res@cnFillOn = True
plot = gsn_csm_contour_map_polar(wks,t1,res)

Being a FORTRAN first user, it takes some getting used to variables
indices starting at zero. There are inconsistencies, because then when
you have functions like the get field the fields start at 1 not zero.
Why not be consistent? (I guess no language is).


“What is wanted is not the will to believe,
but the wish to find out, which is its exact opposite”
  ... Bertrand Russell

On 06/14/2013 09:32 AM, Mary Haley wrote:
> On Jun 11, 2013, at 2:32 PM, Jeremy Winick wrote:
>> I am trying to learn ncl as a long time matlab user who is retired and
>> cannot afford new licenses.
>> My first problem I am attempting is plotting contours of satellite
>> retrieved data on map projections.
>> I have gotten this to work to some extent, but want to automate the
>> process and don't see how to do it efficiently. The data comes from limb
>> retrievals of satellite passes, 15 orbits per day. At least one day is
>> needed to get the longitudinal coverage, and sometimes two or three days
>> will help fill in more points.
>> To contour temperature at some altitude or pressure level, using
>> gsn_csm_contour_map_polar I am currently using natgrid to put the
>> temperature data at each lat, lon point on a uniform grid. Then after
>> setting, units, named dimensions and coordinate variables I can get the
>> plot. With one day of data these contours often have some "non-smooth"
>> features. Is there a way to get smoother contours or better griding.
>> Note that these temperature points are measured at different times and
>> maybe getting a grid should not just use nearest neighbor type methods
>> but is there a way to weight points by how far apart in time?
> Have you tried contouring the data directly, rather than putting it on a grid?
> http://www.ncl.ucar.edu/Applications/contour1d.shtml
> It would help if you could show us some sample images.
>> Now, to the problem of automating. All the examples I see on the web
>> pages surely don't fit my type of need. The data opened in the script
>> uses a specific name. My data is data that has already been analyzed and
>> it in large ascii files - column format. Although all columns are not
>> numeric, readAsciiTable seems to work although the column of the data
>> that are not numeric (say dates in Feb12 type format) are not usable. I
>> want to tell the script what the file to open (changing from run to run
>> or incorporated in an csh-shell or bourne shell script, and also pass
>> which days, and which column (temperature at one of a few levels).
>> Setting enumerable environment variables seems like an awkward solution.
>> Passing these as command line arguments would be very nice (like in
>> matlab), but I haven't seen that this is possible. My work around now
>> would be to open a fixed name file that I read in the script and change
>> it each run to change these important parameters. I think I must be
>> missing something that is basic, since my problem doesn't seem to be
>> unusual.
> It would help if we could see a sample file and script that you've written so far.
> As Jonathan Vigh pointed out, there are different methods for reading ascii
> files, but it just depends on what your file looks like. If you have a mix of data
> and numbers, then I usually find the str_get_field and str_get_cols functions
> to be the best way of reading in columnar data. I use str_get_field if the data
> is separated by some delimiter, like a comma. I use str_get_cols if the data
> are in perfectly lined up columns, which can be specified by column numbers.
> If you can provide your file and script, then you can use our ftp site:
> http://www.ncl.ucar.edu/report_bug.shtml#HowToFTP
>> One other complaint is that the web page support when viewed on Firefox
>> 20.0 does not print correctly. Only the first page (and not even all of
>> it) prints. I know these pages are designed for online work, but having
>> ten or more open at a time and going back to just one or two, it would
>> be easier to have those two printed out.
> Can you give an example of which page is not printing out properly? One issue might be
> that you need to click on the page you want to print, to make sure you are printing the
> text you want, and not the header information.
> --Mary
>> Jeremy
>> --
>> “What is wanted is not the will to believe,
>> but the wish to find out, which is its exact opposite”
>> ... Bertrand Russell
>> _______________________________________________
>> ncl-talk mailing list
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
Received on Fri Jun 14 16:48:21 2013

This archive was generated by hypermail 2.1.8 : Mon Jun 24 2013 - 11:46:47 MDT