[ncl-talk] Problems with latitude dimension
Elaine Silva
laine.usp at gmail.com
Fri Feb 10 19:39:04 MST 2017
Dear ncl users
I'm trying to use the mjoclivar_2.ncl script with daily total precipitation
data from Era Interim, however an error is occurring
Could someone help me? Thank you very much in advance!
Variable: xAnom_sm
Type: float
Total Size: 951857280 bytes
237964320 values
Number of Dimensions: 3
Dimensions and sizes: [time | 1826] x [latitude | 181] x [longitude |
720]
Coordinates:
time: [920436..964236]
latitude: [-45..45]
longitude: [-180..179.5]
Number Of Attributes: 3
_FillValue : -32767
long_name : Anomalies from Smooth Daily Climatology
units : m
(0)
(0) Anomalies from Smooth Daily Climatology: min=-0.0362034
max=0.329606
FileAddVar, in file: NclFile.c, line: 412
fatal:FileAddVar: Dimension (latitude) is not currently defined, can't add
variable
fatal:["Execute.c":8578]:Execute: Error occurred at or near line 130 in
file 01mjoelaine_2_pr.ncl
And the script
; ncl 01mjoelaine_2_pr.ncl
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"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/diagnostics_cam.ncl"
begin
;-------------------------------------------------------------
; User specifications
;-------------------------------------------------------------
NC = True ; Deseja criar o arquivo
netCDF? True = sim e False = não
PLOT = False ; sample plots? Não cria
figura
ymdStrt = 20050101 ; start yyyymmdd
ymdLast = 20091231 ; last
yrStrt = ymdStrt/10000
yrLast = ymdLast/10000
nhar = 4 ; number of fourier comp
var = "tp" ; name of file. Nome da
variável do seu arquivo de entrada
vName = "TP" ; name for plots. Nome da
variável do arquivo netCDF que será gerado
;diri = "/project/cas/shea/CDC/" ; input directory
diri = "../../dados/diario/" ; input directory.
Diretório de entrada dos dados
fili = var+".era.interim.2005.2009.invert.nc" ; input file.
Nome do arquivo a ser lido.
if (NC) then
diro= "./" ; output dir. Onde será
gerado o arquivo netCDF. "./" quer dizer no mesmo diretório do script
;filo= var+".day.anomalies.nc" ; output file
filo= var+".day.anomalies."+yrStrt+"-"+yrLast+".nc" ; output file.
Nome do arquivo netCDF a ser criado.
end if
if (PLOT) then
wksType = "pdf" ; extensão da figura. eps,
pdf, ps, png
wksName = "mjoclivar" ; "mjo."+yrStrt+"_"+yrLast.
Nome da figura a ser gerada.
end if
;***********************************************************
; Read user specified time and create required yyyyddd
;***********************************************************
f = addfile (diri+fili , "r")
time = f->time ; all times on file
ymd = cd_calendar(time, -2) ; yyyymmdd
iStrt = ind(ymd.eq.ymdStrt) ; index start
iLast = ind(ymd.eq.ymdLast) ; index last
delete(time)
delete(ymd)
;***********************************************************
; Read user specified time and create required yyyyddd
;***********************************************************
time = f->time(iStrt:iLast) ; time:units = "hours since"
TIME = cd_calendar(time, 0) ; type float
year = floattointeger( TIME(:,0) )
month = floattointeger( TIME(:,1) )
day = floattointeger( TIME(:,2) )
ddd = day_of_year(year, month, day)
yyyyddd = year*1000 + ddd ; needed for input
;***********************************************************
; Read data: short2flt
;***********************************************************
x = short2flt( f->$var$(iStrt:iLast,:,:) ) ; convert to float
printVarSummary( x )
;***********************************************************
; Compute daily climatology: raw and then 'smoothed'
;***********************************************************
xClmDay = clmDayTLL(x, yyyyddd) ; daily climatology at each grid
point
printVarSummary(xClmDay)
;***********************************************************
; Compute smoothed daily climatology using 'nhar' harmonics
;***********************************************************
xClmDay_sm = smthClmDayTLL(xClmDay, nhar)
printVarSummary(xClmDay_sm)
;***********************************************************
; Compute daily anomalies using raw and smoothed climatologies
;***********************************************************
xAnom = calcDayAnomTLL (x, yyyyddd, xClmDay)
printVarSummary(xAnom)
printMinMax(xAnom, True)
xAnom_sm = calcDayAnomTLL (x, yyyyddd, xClmDay_sm)
xAnom_sm at long_name = "Anomalies from Smooth Daily Climatology"
printVarSummary(xAnom_sm)
printMinMax(xAnom_sm, True)
delete( x ) ; no longer needed
;***********************************************************
; Create netCDF: convenience use 'simple' method
;***********************************************************
dimx = dimsizes(xAnom)
ntim = dimx(0)
nlat = dimx(1)
mlon = dimx(2)
if (NC) then
lat = f->latitude
lon = f->longitude
system("/bin/rm -f "+diro+filo) ; rm any pre-exist file, if any
fnc = addfile (diro+filo, "c")
filAtt = 0
filAtt at title = vName+": Daily Anomalies:
"+yrStrt+"-"+yrLast
filAtt at source_file = fili
filAtt at creation_date = systemfunc("date")
fileattdef( fnc, filAtt ) ; copy file attributes
setfileoption(fnc,"DefineMode",True)
varNC = vName+"_anom"
varNC_sm = vName+"_anom_sm"
dimNames = (/"time", "lat", "lon"/)
dimSizes = (/ -1 , nlat, mlon/)
dimUnlim = (/ True , False, False/)
filedimdef(fnc,dimNames,dimSizes,dimUnlim)
filevardef(fnc, "time" ,typeof(time),getvardims(time))
filevardef(fnc, "lat" ,typeof(lat) ,getvardims(lat))
filevardef(fnc, "lon" ,typeof(lon) ,getvardims(lon))
filevardef(fnc, varNC ,typeof(xAnom) ,getvardims(xAnom))
filevardef(fnc, varNC_sm,typeof(xAnom) ,getvardims(xAnom))
filevarattdef(fnc,"time" ,time) ; copy time
attributes
filevarattdef(fnc,"lat" ,lat) ; copy lat
attributes
filevarattdef(fnc,"lon" ,lon) ; copy lon
attributes
filevarattdef(fnc,varNC ,xAnom)
filevarattdef(fnc,varNC_sm,xAnom_sm)
fnc->time = (/time/)
fnc->lat = (/lat/)
fnc->lon = (/lon/)
fnc->$varNC$ = (/xAnom /)
fnc->$varNC_sm$ = (/xAnom_sm/)
end if
;************************************************
; plotting parameters
;************************************************
if (PLOT) then
LAT = (/ 60, 45, 5, -5, -45, 60 /)
LON = (/270, 30, 90, 90, 180, 0 /)
nPts = dimsizes( LAT )
plot = new ( nPts, graphic)
data = new ( (/2,366/), typeof(xClmDay), getFillValue(xClmDay))
wks = gsn_open_wks (wksType,wksName)
res = True ; plot mods
desired
res at gsnDraw = False
res at gsnFrame = False
res at trXMinF = 1
res at trXMaxF = 366
;res at tiMainString = ""
res at xyLineThicknesses = (/1.0, 2.0/) ; make 2nd lines
thicker
res at xyLineColors = (/"blue","red"/) ; change line
color
res at xyMonoDashPattern = True ; all solid
do np=0,nPts-1
data(0,:) = xClmDay(:,{LAT(np)},{LON(np)})
data(1,:) = xClmDay_sm(:,{LAT(np)},{LON(np)})
res at gsnCenterString = "lat="+LAT(np)+" lon="+LON(np)
plot(np) = gsn_csm_y (wks,data,res) ; create plot
end do
resP = True ; modify the panel
plot
resP at txString = vName+": Raw/Smooth Daily Climatology:
"+yrStrt+"-"+yrLast
resP at gsnMaximize = True
resP at gsnPaperOrientation = "portrait"
gsn_panel(wks,plot,(/(nPts/2),2/),resP) ; now draw as
one plot
;==========
; Plot anomalies for an arbitrarily selected near equatorial location
; Time: Oct 1, 1996 to April 1,1997 [arbitrary selection]
;==========
LATX = 0
LONX = 90
yyyymmdd = cd_calendar(time, -2)
;;yrfrac = yyyymmdd_to_yyyyfrac (yyyymmdd, 0)
;;delete(yrfrac at long_name)
xAnom at long_name = "Anomalies from Raw" ; short labels for plot
xAnom_sm at long_name = "Anomalies from Smooth"
ntBegin = ind(yyyymmdd.eq.20051001)
ntEnd = ind(yyyymmdd.eq.20060401)
monthLabels = (/1,4,7,10/)
monNam = (/"Jan","Feb","Mar","Apr","May","Jun" \
,"Jul","Aug","Sep","Oct","Nov","Dec" /)
xVal = new(ntim, typeof(xAnom&time) , "No_FillValue") ; bigger
than
xLab = new(ntim, "string", "No_FillValue") ; needed
xValm = new(ntim, typeof(xAnom&time) , "No_FillValue") ; bigger
than
ntm = -1
cr = inttochar(10) ; carriage return
do nt=ntBegin,ntEnd
if (day(nt).eq.1) then
ntm = ntm + 1
xVal(ntm) = xAnom&time(nt)
xLab(ntm) = monNam(month(nt)-1)
if (month(nt).eq.1) then
xLab(ntm) = xLab(ntm) + cr +sprinti("%0.4i", year(nt))
end if
end if
end do
rxy = True
rxy at gsnDraw = False
rxy at gsnFrame = False
rxy at gsnYRefLine = 0.0 ; create a reference
line
rxy at gsnAboveYRefLineColor = "red" ; above ref line fill red
rxy at gsnBelowYRefLineColor = "blue" ; below ref line fill blue
rxy at xyLineThicknessF =
2.0
rxy at vpHeightF = 0.4 ; resize
rxy at vpWidthF = 0.8
rxy at tmXBMode = "Explicit"
rxy at tmXBValues = xVal(0:ntm)
rxy at tmXBLabels = xLab(0:ntm)
plot(0) = gsn_csm_xy (wks,time(ntBegin:ntEnd) \
,xAnom(ntBegin:ntEnd,{0},{90}),rxy)
plot(1) = gsn_csm_xy (wks,time(ntBegin:ntEnd) \
,xAnom_sm(ntBegin:ntEnd,{0},{90}),rxy)
resP at txString = vName+": Anomalies: (0,90E)"
gsn_panel(wks,plot(0:1),(/2,1/),resP)
end if
end
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20170211/59a70fb1/attachment.html
More information about the ncl-talk
mailing list