[ncl-talk] Hovmoller 2d lat and 2d lon rectangular grid
David Adams
dave.k.adams at gmail.com
Tue Sep 15 15:37:42 MDT 2015
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
>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150915/19cf7baa/attachment.html
More information about the ncl-talk
mailing list