[ncl-talk] Problems with latitude dimension
Dennis Shea
shea at ucar.edu
Sat Feb 11 08:03:26 MST 2017
Possibly .....
After
lat = f->latitude
lon = f->longitude
Add
lat!0 = "lat"
lon!0 = "lon"
Good Luck
On Fri, Feb 10, 2017 at 7:39 PM, Elaine Silva <laine.usp at gmail.com> wrote:
> 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
>
>
>
> _______________________________________________
> 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/20170211/0d8e4eb6/attachment.html
More information about the ncl-talk
mailing list