[ncl-talk] Hovmoller 2d lat and 2d lon rectangular grid
Dennis Shea
shea at ucar.edu
Wed Sep 16 08:59:47 MDT 2015
General comment:
You can do a 'cross-section' on a curvilinear grid. However, it
depends on what you mean by cross-section ... :-)
Rectilinear grid: 'classic' hovmueller-sections: x(time,lat, lon)
(a) x(:,nl,:) ; west-east at the specified latitude index
=> x(time,longitude)
(b) x(:,:,ml) ; south-north (or vice-versa) at a specified
longitude index => x(time,latitude)
In case (a) the only values only are needed for the (say)
x-axis; case (b) only the lat values
Curvilinear grid: y(time,nlat,mlon) ... generically the same but you
must plot lat/lon pair. Although not Hovmueller plots. See:
http://www.ncl.ucar.edu/Applications/narr.shtml
Example 6: These show pressure levels on the left but they could be time
On Tue, Sep 15, 2015 at 3:37 PM, David Adams <dave.k.adams at gmail.com> wrote:
> Hi Mary et al.
> sorry for the delay. I have got the Hovmoller code to work for the GOES IR
> Brightness temperature data. However, what surprises me
> is that it works with the original curvilinear data, which it wasn´t
> supposed to. That is, this GOES data is not on a rectangular grid.
>
> Here is the code
>
> iri = "./"
> fili = systemfunc("ls goes13_4_2013_190*subset.nc")
> ; print(fili)
> nfili= dimsizes(fili)
> ; print("nfili="+nfili)
>
> do i=0,nfili-1
> f = addfile(fili(i),"r")
> T = f->new_temp ; float data(time, yc, xc) ;
> data:type = "VISR" ;
> la = f->new_lat
> lo = f->new_lon
> nT = dimsizes(T)
> if (i .eq. 0)
> d=new((/nfili, nT(0), nT(1)/), float)
> lat=d
> lon=d
> d at _FillValue = 0 ; not on file
> d!0="Time"
> end if
> d(i,:,:) = T
> lat(i,:,:) = la
> lon(i,:,:) = lo
> delete([/T, la, lo/])
> end do
>
> ; printVarSummary(d)
> ; print("d: min="+min(d)+" max="+max(d))
>
> nl = 180
> dhov = d(:,nl,:) ; (time,longitude) ; <== plot
> print(dhov)
> d at lat2d = lat(:,nl,:) ; (yc, xc)
> d at lon2d = lon(:,nl,:)
>
> ;*********************************
> ; create plot
> ;*********************************
> fNameBase = str_get_field(fili,1,".") + str_get_field(fili,2,".")
>
> pltName = "goes_13"
> pltType = "pdf" ; "ps", "eps", "pdf",
> "png","x11"
> pltDir = "./"
>
> wks = gsn_open_wks(pltType, pltDir+pltName)
> gsn_define_colormap(wks,"BlAqGrYeOrReVi200"); specify a color map
>
> res = True
> res at cnFillOn = True ; turn on color
> res at cnFillMode = "RasterFill" ; cell mode
> res at cnLinesOn = False ; Turn off contour lines
>
> res at sfXArray = lon(0,nl,:)
>
> res at gsnSpreadColors = True ; use full colormap
> res at gsnAddCyclic = False ; data not cyclic
> res at gsnMaximize = True ; ps, pdf, pdf
> res at pmTickMarkDisplayMode = "Always" ; use NCL default
> res at lbLabelAutoStride = True ; let NCL decide spacing
> res at cnLevelSelectionMode = "ManualLevels" ;"ExplicitLevels"
> res at cnMinLevelValF = 180.0 ; set the
> minimum contour level
> res at cnMaxLevelValF = 310. ; set the
> maximum contour level
> res at cnLevelSpacingF = 5.0 ; set the
> contour interval
> res at cnRasterSmoothingOn = True
> res at lbLabelStride = 5.0 ; every other label bar label
>
> plot = gsn_csm_contour(wks,dhov, res)
>
> -------------------------PrintSummary-------------------------------------
> Variable: dhov
> Type: float
> Total Size: 48504 bytes
> 12126 values
> Number of Dimensions: 2
> Dimensions and sizes: [Time | 86] x [xc | 141]
> Coordinates:
> Number Of Attributes: 6
> _FillValue : 1
> time : 1376178300
> units : K
> coordinates : lon lat
> type : GVAR
> long_name : Temperature
>
> thanks,
> Dave
>
>
>
>
>
>
> On Mon, Aug 17, 2015 at 11:44 PM, Mary Haley <haley at ucar.edu> wrote:
>>
>> Dave,
>>
>> You are plotting "dhov", which according to your comments is time x
>> longitude. The error is:
>>
>> warning:IrTransInitialize: error creating spline approximation for
>> trYCoordPoints; defaulting to linear
>>
>> which is referencing the Y coordinate array, and hence the time array, in
>> your case.
>>
>> Usually if there's a spline error, it's because the coordinate points are
>> too irregular. That is, they might be closely spaced in one section of your
>> data, and then widely spaced in others. This is very true of time arrays,
>> if your time array is not constructed properly.
>>
>> This is a case where you need to look at your data, and in particular, the
>> "time" array that's attached to dhov. What does the following report:
>>
>> printVarSummary(dhov)
>> print(dhov&time)
>>
>> --Mary
>>
>>
>> On Mon, Aug 17, 2015 at 10:51 AM, David Adams <dave.k.adams at gmail.com>
>> wrote:
>>>
>>> Hi NCLers,
>>>
>>> I have managed to get GOES IR data from curvilinear to rectangular grids.
>>> Lat and Lon are still 2d, however. But the code I have for the Hovmoller
>>> diagrams gives this error. I am sure I am just not reading the data
>>> correctly.
>>>
>>> The error, ncl code and ncdump of the file follow:
>>>
>>> thanks,
>>> Dave
>>>
>>> ---------------------- error
>>> ------------------------------------------------------
>>> warning:IrTransInitialize: error creating spline approximation for
>>> trYCoordPoints; defaulting to linear
>>> fatal:ContourPlotDraw: CPCICA - ONE OF THE CORNER POINTS OF THE CELL
>>> ARRAY IS INCORRECT
>>> fatal:ContourPlotDraw: draw error
>>> fatal:ContourPlotDraw: draw error
>>> fatal:PlotManagerDraw: error in plot draw
>>> fatal:_NhlPlotManagerDraw: Draw error
>>>
>>> ------------------------ NCL Code Snippet
>>> -----------------------------------------------
>>> diri = "./"
>>> fili = systemfunc("ls goes1*.nc")
>>>
>>> print(fili)
>>> nfili= dimsizes(fili)
>>> print("nfili="+nfili)
>>>
>>> f = addfiles(diri+fili,"r")
>>> d = f[:]->data ; float data(time, yc, xc) ; data:type
>>> = "IR" ;
>>> d at _FillValue = 330
>>> lat = f[:]->lat
>>> lon = f[:]->lon
>>> time= f[:]->time
>>>
>>> nl = 600
>>> dhov = d(:,nl,:) ; (time,longitude) ; <== plot
>>> ; Fix the variable so it has recognizable missing data
>>> ;*********************************
>>>
>>> lat2d = f[:]->lat
>>> lon2d = f[:]->lon
>>> d at lat2d = lat2d ; (yc, xc)
>>> d at lon2d = lon2d
>>>
>>> crTime = f[:]->crTime ; assorted info
>>> crDate = f[:]->crDate
>>> lineRes = f[:]->lineRes
>>> elemRes = f[:]->elemRes
>>>
>>> ;*********************************
>>> ; create plot
>>> ;*********************************
>>> fNameBase = str_get_field(fili,1,".") + str_get_field(fili,2,".")
>>>
>>> ;pltName = fNameBase
>>> pltName = "goes_4"
>>> pltType = "pdf" ; "ps", "eps", "pdf",
>>> "png","x11"
>>> pltDir = "./"
>>>
>>> wks = gsn_open_wks(pltType, pltDir+pltName)
>>> gsn_define_colormap(wks,"BlAqGrYeOrReVi200"); specify a color map
>>>
>>> res = True
>>> res at cnFillOn = True ; turn on color
>>> res at cnFillMode = "RasterFill" ; cell mode
>>> res at cnLinesOn = False ; Turn off contour line
>>> ; Create longitude labels on the x-axis
>>> res at gsnSpreadColors = True ; use full colormap
>>> res at gsnAddCyclic = False ; data not cyclic
>>> res at gsnMaximize = True ; ps, pdf, pdf
>>> res at pmTickMarkDisplayMode = "Always" ; use NCL default
>>> res at lbLabelAutoStride = True ; let NCL decide spacing
>>> res at cnLevelSelectionMode = "ManualLevels" ;"ExplicitLevels"
>>> res at cnMinLevelValF = 180.0 ; set the
>>> minimum contour level
>>> res at cnMaxLevelValF = 310. ; set the
>>> maximum contour level
>>> res at cnLevelSpacingF = 5.0 ; set the
>>> contour interval
>>> res at cnRasterSmoothingOn = True
>>> res at lbLabelStride = 5.0 ; every other label bar
>>> label
>>>
>>> plot = gsn_csm_contour(wks,dhov(:,:), res)
>>>
>>> ---------------------------------------------------------------------------------------------------------
>>>
>>> netcdf goes13_4_2013_210_2002_1byte_rect {
>>> dimensions:
>>> xc = 2720 ;
>>> yc = 1301 ;
>>> time = 1 ;
>>> auditCount = 2 ;
>>> auditSize = 80 ;
>>> variables:
>>> int version ;
>>> version:long_name = "McIDAS area file version" ;
>>> int sensorID ;
>>> sensorID:long_name = "McIDAS sensor number" ;
>>> int imageDate ;
>>> imageDate:long_name = "image year and day of year (in ccyyddd
>>> format)" ;
>>> int imageTime ;
>>> imageTime:long_name = "image time in UTC (hours/minutes/seconds,
>>> in HHMMSS format)" ;
>>> int startLine ;
>>> startLine:long_name = "starting image line (in satellite
>>> coordinates)" ;
>>> int startElem ;
>>> startElem:long_name = "starting image element (in satellite
>>> coordinates)" ;
>>> int time(time) ;
>>> time:long_name = "seconds since 1970-1-1 0:0:0" ;
>>> time:units = "seconds since 1970-1-1 0:0:0" ;
>>> int dataWidth ;
>>> dataWidth:long_name = "number of 8-bit bytes per source data
>>> point" ;
>>> int lineRes ;
>>> lineRes:long_name = "resolution of each pixel in line direction"
>>> ;
>>> lineRes:units = "km" ;
>>> int elemRes ;
>>> elemRes:long_name = "resolution of each pixel in element
>>> direction" ;
>>> elemRes:units = "km" ;
>>> int prefixSize ;
>>> prefixSize:long_name = "line prefix size in 8-bit bytes" ;
>>> int crDate ;
>>> crDate:long_name = "image creation year and day of year in
>>> ccyyddd format" ;
>>> int crTime ;
>>> crTime:long_name = "image creation time in UTC in hhmmss format"
>>> ;
>>> int bands ;
>>> bands:long_name = "satellite channel number" ;
>>> char auditTrail(auditCount, auditSize) ;
>>> auditTrail:long_name = "audit trail" ;
>>> float data(time, yc, xc) ;
>>> data:long_name = "Temperature" ;
>>> data:type = "VISR" ;
>>> data:coordinates = "lon lat" ;
>>> data:units = "K" ;
>>> float lat(yc, xc) ;
>>> lat:long_name = "lat" ;
>>> lat:units = "degrees_north" ;
>>> float lon(yc, xc) ;
>>> lon:long_name = "lon" ;
>>> lon:units = "degrees_east" ;
>>>
>>> // global attributes:
>>> :Conventions = "CF-1.4" ;
>>> :Source = "McIDAS Area File" ;
>>> :Satellite\ Sensor = "G-13 IMG " ;
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> ncl-talk at ucar.edu
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
More information about the ncl-talk
mailing list