Helen,
The error messages are telling you quite a bit here. For one:
warning:Argument 4 of the current function or procedure was coerced to the
appropriate type and thus will not change if the function or procedure
Your variables are probably coming in as floats, and you are trying to
pass these floats to a Fortran routine that's expecting doubles.
You should convert your floats to doubles if you want to get rid
of these warnings.
Secondly, the error:
fatal:VENUSDRV: dimension size of dimension (3) of V must be equal to the
value of NTIM
indicates that your arrays do not have the correct dimensions. I
believe this is b/c you are passing the arrays in the wrong order.
Fortran and NCL order their arrays differently (column order versus
row order). This means if you have an array on the Fortran
side that is ntime x nlev x nlat x nlon, then in the NCL
script, that array must be ordered nlon x nlat x nlev x ntime.
In your case, the NCL's "v" is
time x lev_p lat x lon
and the Fortran routine is expecting
V(NTIM,KLEV,NLAT,MLON)
You will either need to reorder the arrays on the Fortran side,
or on the NCL side:
vreord = v(lon|:,lat|:,lev_P|:,time|:)
You have to reorder "meridstream" and probably PS as well:
PSreord = PS(lon|:,lat|:,time|:)
meridstream = new ( (/nlat,klvl,ntim/) , typeof(vavg), getFillValue(vavg))
HELEN::VENUSDRV (mlon,nlat,klvl,ntim,vreord,lat,lev_p,PSreord, \
vavg@_FillValue,meridstream)
You could instead reorder the arrays on the Fortran side. This will
probably run faster than reordering them on the NCL side.
--Mary
On Thu, 2 Jul 2009, Helen Parish wrote:
> I am trying to plot the stream function for my data, using an NCL script. I
> am getting the error below, but I do not see what is causing it and how to
> get the NCL script to run.
>
> I have included the NCL script that produces the error and the Fortran
> routine that the script calls (which had to be customized for my model from
> DZPSIDRV). I also include the statement I use to compile the Fortran code.
>
> Does anyone know what is causing this error and how to get the script to run
> ?.
>
> Thanks,
> Helen.
>
>
> a) The error:
>
> Variable: v
> Type: float
> Total Size: 344448 bytes
> 86112 values
> Number of Dimensions: 4
> Dimensions and sizes: [time | 1] x [lev_p | 26] x [lat | 46] x [lon | 72]
> Coordinates:
> time: [13710..13710]
> lev_p: [91500..0.03]
> lat: [ -90.. 90]
> lon: [ 0.. 355]
> Number Of Attributes: 4
> units : m/s
> long_name : Meridional wind
> cell_method : time: mean
> _FillValue : -9999
> warning:Argument 4 of the current function or procedure was coerced to the
> appropriate type and thus will not change if the function or procedure
> modifies its value
> warning:Argument 6 of the current function or procedure was coerced to the
> appropriate type and thus will not change if the function or procedure
> modifies its value
> warning:Argument 7 of the current function or procedure was coerced to the
> appropriate type and thus will not change if the function or procedure
> modifies its value
> warning:Argument 9 of the current function or procedure was coerced to the
> appropriate type and thus will not change if the function or procedure
> modifies its value
> fatal:VENUSDRV: dimension size of dimension (3) of V must be equal to the
> value of NTIM
> fatal:Execute: Error occurred at or near line 66 in file zonstreamtest.ncl
>
>
> b) The NCL script:
>
> external HELEN "./venus_mpsitest.so"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
> begin
> ; input file
> diri = "./"
> fili = "surfha.cam2.h0.0252-08-25-00000.nc"
> fi = addfile (diri+fili, "r")
> ; output file
> diro = "./"
> filo = "helen.nc"
> system ("/bin/rm -f "+diro+filo) ; remove any pre-exist file
> fo = addfile (diro+filo, "c")
> ; add any file attributes
> fo_at_title = fi_at_title + ": Selected Variables at pressure levels"
> fo_at_history = systemfunc ("date")
> fo_at_source = fi_at_source
> fo_at_case = fi_at_case
> fo_at_Conventions = fi_at_Conventions
>
> Var = (/ "T" , "Q", "U", "V"/) ; select variable to be interpolated.
> nVar = dimsizes (Var)
> ; desired output levels
> lev_p = (/ 91500., 90000., 87500., 80000., 67500., 50000., 30000., 20000.,
> 10000., 5000., 3000., 1500., 750., 390., 200., 100., 50.,25., 14., 7., 3.5,
> 1.9, 0.95, 0.4, 0.1, 0.03 /)
> lev_p!0 = "lev_p" ; variable and dimension name the same
> lev_p&lev_p = lev_p ; create coordinate variable
> lev_p_at_long_name = "pressure" ; attach some attributes
> lev_p_at_units = "hPa"
> lev_p_at_positive = "down"
>
> hyam = fi->hyam ; read hybrid info
> hybm = fi->hybm
> PS = fi->PS
> P0mb = 0.01*fi->P0
>
> do n=0,nVar-1 ; loop over the variables
> X = fi->$Var(n)$
> Xp = vinth2p (X, hyam, hybm, lev_p ,PS, 1, P0mb, 2, True)
> copy_VarAtts(X, Xp)
> fo->$Var(n)$ = Xp ; write to netCDF file
> print (Var(n)+": interpolated and written to netCDF")
> end do
>
> end
> ;***********************
> ; zonal.ncl
> ;***********************
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> ;***********************
> begin
> f = addfile ("helen.nc","r")
> print(f)
> v = f->V
> lat = f->lat
> printVarSummary( v )
>
> dimu = dimsizes( v )
> ntim = dimu(0)
> klvl = dimu(1)
> nlat = dimu(2)
> mlon = dimu(3)
>
> vavg = dim_avg_Wrap(v)
>
> meridstream = new ( (/ntim,klvl,nlat/) , typeof(vavg), getFillValue(vavg))
> HELEN::VENUSDRV
> (mlon,nlat,klvl,ntim,v,lat,lev_p,PS,vavg@_FillValue,meridstream)
>
> ;***********************
> ; Create Plot
> ;***********************
> wks = gsn_open_wks ("pdf", "surfha2520825stream" ) ; open
> workstation
> gsn_define_colormap(wks,"rainbow") ; choose colormap
>
> res = True
> res_at_cnFillOn = True
> res_at_lbLabelAutoStride = True
> res_at_gsnMaximize = True ; if [ps, eps, pdf] make large
> res_at_gsnSpreadColors = True ; span color map
>
> do nt=0,ntim-1
>
> res_at_gsnCenterString = "surfha2520825stream"
> plot = gsn_csm_pres_hgt(wks, meridstream(nt,:,:), res ) ; (lev,lat)
>
> res_at_trYReverse = True
>
> kl = 5
> res_at_mpFillOn = False
> res_at_mpGridAndLimbOn = True
> res_at_mpGridLineDashPattern = 2
> res_at_mpOutlineBoundarySets = "NoBoundaries"
> res_at_mpCenterLonF = 180.
> res_at_gsnCenterString = res_at_gsnCenterString+" p="+u&lev_p(kl)
>
> end do
>
> end
>
> c) The Fortran routine VENUSDRV :
>
> C NCLFORTSTART
> SUBROUTINE VENUSDRV(MLON,NLAT,KLEV,NTIM,V,LAT,P,PS,VMSG,
> + ZMPSI)
> IMPLICIT NONE
>
> c NCL: zmpsi = zon_mpsi(v,lat,p,ps)
> c INPUT
> INTEGER MLON,NLAT,KLEV,NTIM
>
> DOUBLE PRECISION V(NTIM,KLEV,NLAT,MLON),LAT(NLAT),P(KLEV),
> + PS(NTIM,NLAT,MLON),VMSG
> c OUTPUT
> DOUBLE PRECISION ZMPSI(NTIM,KLEV,NLAT)
> C NCLEND
> c LOCAL
> DOUBLE PRECISION PTMP(2*KLEV+1),DP(2*KLEV+1),VVPROF(2*KLEV+1)
> DOUBLE PRECISION VBAR(NTIM,KLEV,NLAT),VTMP(NTIM,KLEV,NLAT,MLON)
>
> CALL VENUSZON(MLON,NLAT,KLEV,NTIM,V,LAT,P,PS,VMSG,ZMPSI,
> + PTMP,DP,VTMP,VBAR,VVPROF)
>
> RETURN
> END
> c --
> SUBROUTINE VENUSZON(MLON,NLAT,KLEV,NTIM,V,LAT,P,PS,VMSG,ZMPSI,
> + PTMP,DP,VTMP,VBAR,VVPROF)
> IMPLICIT NONE
>
> INTEGER MLON,NLAT,KLEV,NTIM
> c DOUBLE PRECISION V(MLON,NLAT,KLEV,NTIM),LAT(NLAT),P(KLEV),
> c + PS(MLON,NLAT,NTIM),VMSG
> DOUBLE PRECISION V(NTIM,KLEV,NLAT,MLON),LAT(NLAT),P(KLEV),
> + PS(NTIM,NLAT,MLON),VMSG
> c output
> c DOUBLE PRECISION ZMPSI(NLAT,KLEV)
> c DOUBLE PRECISION ZMPSI(KLEV,NLAT)
> DOUBLE PRECISION ZMPSI(NTIM,KLEV,NLAT)
> c temporary / local
> DOUBLE PRECISION PTMP(0:2*KLEV),DP(0:2*KLEV),VVPROF(0:2*KLEV)
> DOUBLE PRECISION VTMP(NTIM,KLEV,NLAT,MLON),VBAR(NTIM,KLEV,NLAT)
> c DOUBLE PRECISION VTMP(MLON,NLAT,KLEV,NTIM),VBAR(NLAT,KLEV)
>
> DOUBLE PRECISION PI,G,A,C,RAD,CON,PTOP,PBOT
> INTEGER ML,NL,KL,NT,KNT,KFLAG
>
> G = 9.80616D0
> A = 6.37122D6
> PI = 4.D0*ATAN(1.D0)
> RAD = PI/180.D0
> CON = 2.D0*PI*A/G
>
> PTOP = 3.D0
> PBOT = 9150000.D0
>
> c calculate press at all levels [even half levels]
>
> KNT = 0
> DO KL = 1,2*KLEV - 1,2
> KNT = KNT + 1
> PTMP(KL) = P(KNT)
> END DO
>
> DO KL = 2,2*KLEV - 2,2
> PTMP(KL) = (PTMP(KL+1)+PTMP(KL-1))*0.5D0
> END DO
>
> PTMP(0) = PTOP
> PTMP(2*KLEV) = PBOT
>
> c dp at all levels
>
> DP(0) = 0.D0
> DO KL = 1,2*KLEV - 1
> DP(KL) = PTMP(KL+1) - PTMP(KL-1)
> END DO
> DP(2*KLEV) = 0.D0
>
> c make copy; set any p > ps to msg
>
> KNT = 0
> DO ML = 1,MLON
> DO NL = 1,NLAT
> DO KL = 1,KLEV
> DO NT = 1,NTIM
> IF (P(KL).GT.PS(NT,NL,ML)) THEN
> KNT = KNT + 1
> VTMP(NT,KL,NL,ML) = VMSG
> ELSE
> VTMP(NT,KL,NL,ML) = V(NT,KL,NL,ML)
> END IF
> END DO
> END DO
> END DO
> END DO
>
> c compute zonal mean v using the vtmp variable
>
> DO NT = 1,NTIM
> DO NL = 1,NLAT
> DO KL = 1,KLEV
> KNT = 0
> VBAR(NT,KL,NL) = 0.0D0
> DO ML = 1,MLON
> IF (VTMP(NT,KL,NL,ML).NE.VMSG) THEN
> KNT = KNT + 1
> VBAR(NT,KL,NL) = VBAR(NT,KL,NL) + VTMP(NT,KL,NL,ML)
> END IF
> END DO
> IF (KNT.GT.0) THEN
> VBAR(NT,KL,NL) = VBAR(NT,KL,NL)/DBLE(KNT)
> ELSE
> VBAR(NT,KL,NL) = VMSG
> END IF
> END DO
> END DO
> END DO
>
> c compute mpsi at each latitude [reuse ptmp]
>
> DO NT = 1,NTIM
> DO NL = 1,NLAT
> C = CON*COS(LAT(NL)*RAD)
>
> PTMP(0) = 0.0D0
> DO KL = 1,2*KLEV
> PTMP(KL) = VMSG
> END DO
> DO KL = 0,2*KLEV,2
> VVPROF(KL) = 0.0D0
> END DO
>
> KNT = 0
> DO KL = 1,2*KLEV - 1,2
> KNT = KNT + 1
> VVPROF(KL) = VBAR(NT,KNT,NL)
> END DO
>
> c integrate from top of atmosphere down for each level where vbar
> c is not missing
>
> DO KL = 1,2*KLEV - 1,2
> KFLAG = KL
> IF (VVPROF(KL).EQ.VMSG) GO TO 10
> PTMP(KL+1) = PTMP(KL-1) - C*VVPROF(KL)*DP(KL)
> END DO
> 10 CONTINUE
>
> c impose lower boundary condition to ensure the zmpsi is 0
> c at the bottom boundary
>
> PTMP(KFLAG+1) = -PTMP(KFLAG-1)
>
> c streamfunction is obtained as average from its' values
> c at the intermediate half levels. The minus sign before
> c PTMP in the last loop is to conform to a CSM convention.
>
> DO KL = 1,KFLAG,2
> PTMP(KL) = (PTMP(KL+1)+PTMP(KL-1))*0.5D0
> END DO
>
> KNT = 0
> DO KL = 1,2*KLEV - 1,2
> KNT = KNT + 1
> IF(PTMP(KL).NE.VMSG) THEN
> ZMPSI(NT,KNT,NL) = -PTMP(KL)
> ELSE
> ZMPSI(NT,KNT,NL) = PTMP(KL)
> END IF
> END DO
> END DO
> END DO
>
> RETURN
> END
>
>
> d) Statement used to compile the Fortran code :
>
> WRAPIT -g95 -L /Users/helenparish/g95-install/lib/gcc-lib/powerpc-
> apple-darwin6.8/4.0.3 -l f95 venus_mpsitest.f
>
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Jul 02 2009 - 11:25:48 MDT
This archive was generated by hypermail 2.2.0 : Thu Jul 02 2009 - 11:39:38 MDT