[ncl-talk] multiple ascii file outputs

Mary Haley haley at ucar.edu
Mon May 9 10:40:22 MDT 2016


Please respond back to ncl-talk if you have follow-up questions.

You have two sets of double do loops, where the first set of loops reads in
the data, and the second loop tries to write it.

The problem is that all the data you're reading in the first set of loops
gets reset every time, so by the time you get to your second set of do
loops, you don't have any of the original data left to write.

You want to use just one set of do loops, so that you read in your data,
and then immediately write it to a file before you then read the next set
of data.

I moved all your 'wrf_user_getvar' calls to the outside do loop, because
they only change for every value of "it" and not for every value of "ip".
It's time-consuming to have code in a loop that you don't need, so it's
important to keep your "do" loops as clutter free as possible.

This is UNTESTED:


 ; --------------  LOAD FUNCTIONS AND PROCEDURES ----------------

    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/wrf/WRFUserARW.ncl"

    ; --------------  BEGINING OF NCL SCRIPT ----------------
begin
DATADir = "/mnt/gpfs1/data/private/fgvy024/Build_WRF/WRFV6/run/"
FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d03* ")
numFILES = dimsizes(FILES)
a = addfiles(FILES+".nc","r")
times = wrf_user_list_times(a)     ; get times in the file
ntimes = dimsizes(times)           ; number of times in the file or 1584
for domain1 4750

;new parameters naming
temp = new(ntimes,float)
rain_expp = new(ntimes,float)
rain_conn = new(ntimes,float)
rain_tot = new(ntimes,float)
wspd = new(ntimes,float)
wdir = new(ntimes,float)
snoww_eq = new(ntimes,float)
snowh_h = new(ntimes,float)
;rhh = new(ntimes,float)

ip_lats = (/ 27.3081, 27.6331,27.86,27.959167, \
             27.696667, 27.8024,27.8953, 27.99,27.9818, \
             27.90506667, 27.89671667, 27.901, 27.8987, \
             27.89261111, 27.87828333, 27.85911111, \
             27.60111, 27.74661, 27.85722 /)
ip_lons = (/ 86.5042, 86.2331, 86.46, 86.81278, 86.721389, 86.7144, \
             86.8188, 86.83, 86.7649, 86.37536667, 86.37428333, 86.37625,
86.37911667, \
             86.37594444,86.43361667, 86.43388889, 86.74028, 86.713,
86.79417 /)
res = True
res at returnInt = True

data = new(ntimes, "string")

;******************************* time loop
do it = 0, ntimes-1
  tc2      = wrf_user_getvar(a,"T2",it)       ; temp extraction
  tc2      = tc2 -273.16
  rain_exp = wrf_user_getvar(a,"RAINNC",it)
  rain_con = wrf_user_getvar(a,"RAINC",it)
  u10      = wrf_user_getvar(a,"U10",it)    ; u at 10 m, mass point
  v10      = wrf_user_getvar(a,"V10",it)    ; v at 10 m, mass point
  snow_eq  = wrf_user_getvar(a,"SNOW",it)
  snow_h   = wrf_user_getvar(a,"SNOWH",it)

 ;************************************************
  do ip = 0, 18             ;loop for number of station locations
    point = wrf_user_ll_to_ij(a, ip_lons(ip), ip_lats(ip), True)
    x = point(0) - 1
    y = point(1) - 1

    temp(it)      = tc2(y,x) ; temperature extraction for station location
    rain_expp(it) = rain_exp(y,x)
    rain_conn(it) = rain_con(y,x)
    rain_tot(it)  = rain_exp(y,x) + rain_con(y,x)
    wspd(it)      = ndtooned( sqrt(u10(y,x)^2 + v10(y,x)^2) )
    wdir(it)      = ndtooned( atan2(u10(y,x),v10(y,x))/0.01745329 +180. )
    snoww_eq(it)  = snow_eq(y,x)
    snowh_h (it)  = snow_h(y,x)

;---Create filename to write to
    fName =
ip+"T2oben_1kmtemp_rainnonc_rainc_raintot_wspd_wdir_snow_snowh.txt"

    data (it)= (sprintf("%5.0f",it)    +" " \
          +sprintf("%21.2f", temp(it)) +"  " \
          +sprintf("%22.2f", rain_expp(it)) +"  " \
          +sprintf("%23.2f", rain_conn(it)) +"  " \
          +sprintf("%24.2f", rain_tot(it)) +"  " \
          +sprintf("%17.2f", wspd(it)) +"  " \
          +sprintf("%25.2f", wdir(it)) +"  " \
          +sprintf("%18.2f", snoww_eq(it)) +"  " \
          +sprintf("%19.2f", snowh_h(it)) +"  " )
  end do                      ; end of "ip" loop
  asciiwrite (fName , data)
end do                        ; end of time loop
end


On Mon, May 9, 2016 at 10:23 AM, Ramchandra Karki <rammetro at hotmail.com>
wrote:

> Dear Mary
> thank you very much. i am new to these things. the script is not giving
> errors but i am not getting what i am expecting. i want to get
> -individual ascii file outputs for each stations (19 ascii files for 19
> stations containing time series for whole period). would you please kindly
> suggest me where i need to change to get it?
>
>
> *Regards*
> *ram*
>
>
> ------------------------------
> Date: Mon, 9 May 2016 10:13:48 -0600
> Subject: Re: [ncl-talk] multiple ascii file outputs
> From: haley at ucar.edu
> To: rammetro at hotmail.com
> CC: ncl-talk at ucar.edu
>
>
> It helps if you provide more information about what the problem is.  Are
> you getting an error of any kind? If so, please indicate which line the
> error occurs on.
>
> I think the problem is here:
>
>    do  ip = 0, 18
>     do it = 0,ntimes-1 ;ntimes-1
>     data (it)= (sprintf("%5.0f",it)    +" " \
>           +sprintf("%21.2f", temp(it)) +"  " \
>           +sprintf("%22.2f", rain_expp(it)) +"  " \
>           +sprintf("%23.2f", rain_conn(it)) +"  " \
>           +sprintf("%24.2f", rain_tot(it)) +"  " \
>           +sprintf("%17.2f", wspd(it)) +"  " \
>           +sprintf("%25.2f", wdir(it)) +"  " \
>           +sprintf("%18.2f", snoww_eq(it)) +"  " \
>           +sprintf("%19.2f", snowh_h(it)) +"  " )
> end do
>  end do                     ; end of time loop
>
>  asciiwrite (fName , data)
>
> You have a double "do" loop, where "data" is getting created inside the
> inner do loop.  For each iteration of "ip", then, the previous strings
> assigned to "data" are getting clobbered.  In fact, I'm not even sure why
> you are looping across "ip=0,18", because the inside loop is doing the
> exact same thing every time.
>
> Finally, you call "asciiwrite", but this is only going to contain "data"
> from the last iteration of "ip" (ip=18).
> --Mary
>
>
>
> On Sat, May 7, 2016 at 8:52 AM, Ramchandra Karki <rammetro at hotmail.com>
> wrote:
>
>  Dear NCL help,
> I am trying to extract multiple stations time series (ascci output for
> each stations) from wrfoutput file but i could not succeed in doing so.
> Would you please kindly guide me to get the solution?
>
>  ; --------------  LOAD FUNCTIONS AND PROCEDURES ----------------
>
>     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/wrf/WRFUserARW.ncl"
>
>     ; --------------  BEGINING OF NCL SCRIPT ----------------
> begin
> DATADir = "/mnt/gpfs1/data/private/fgvy024/Build_WRF/WRFV6/run/"
> FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d03* ")
> numFILES = dimsizes(FILES)
> a = addfiles(FILES+".nc","r")
> times = wrf_user_list_times(a)     ; get times in the file
> ntimes = dimsizes(times)           ; number of times in the file or 1584
> for domain1 4750
>
> ;new parameters naming
> temp = new(ntimes,float)
> rain_expp = new(ntimes,float)
> rain_conn = new(ntimes,float)
> rain_tot = new(ntimes,float)
> wspd = new(ntimes,float)
> wdir = new(ntimes,float)
> snoww_eq = new(ntimes,float)
> snowh_h = new(ntimes,float)
> ;rhh = new(ntimes,float)
>
> ;******************************* time loop
> do it = 0, ntimes-1
> time = it
>
>  ;************************************************
> do ip = 0, 18             ;loop for number of station locations
> ip_lats = (/ 27.3081, 27.6331, 27.86, 27.959167, \
>           27.696667, 27.8024, 27.8953, 27.99, 27.9818, \
>           27.90506667, 27.89671667, 27.901, 27.8987, \
>           27.89261111, 27.87828333, 27.85911111, \
>           27.60111, 27.74661, 27.85722 /)
> ip_lons = (/ 86.5042, 86.2331, 86.46, 86.81278, 86.721389, 86.7144, \
>            86.8188, 86.83, 86.7649, 86.37536667, 86.37428333, 86.37625,
> 86.37911667, \
>    86.37594444, 86.43361667, 86.43388889, 86.74028, 86.713, 86.79417 /)
> res = True
> res at returnInt = True
> point = wrf_user_ll_to_ij(a, ip_lons(ip), ip_lats(ip), True)
> x = point(0) - 1
> y = point(1) - 1
>
> tc2 = wrf_user_getvar(a,"T2",it)       ; temp extraction
> tc2 = tc2 -273.16
> temp(it) = tc2(y,x) ; temperature extraction for station location
> rain_exp = wrf_user_getvar(a,"RAINNC",it)
> rain_con = wrf_user_getvar(a,"RAINC",it)
> rain_expp(it) = rain_exp(y,x)
> rain_conn(it) = rain_con(y,x)
> rain_tot(it) = rain_exp(y,x) + rain_con(y,x)
> u10 = wrf_user_getvar(a,"U10",it)    ; u at 10 m, mass point
> v10 = wrf_user_getvar(a,"V10",it)    ; v at 10 m, mass point
> wspd(it) = ndtooned( sqrt(u10(y,x)^2 + v10(y,x)^2) )
> wdir(it) = ndtooned( atan2(u10(y,x),v10(y,x))/0.01745329 +180. )
> snow_eq = wrf_user_getvar(a,"SNOW",it)
> snow_h = wrf_user_getvar(a,"SNOWH",it)
> snoww_eq(it) = snow_eq(y,x)
> snowh_h (it) = snow_h(y,x)
>
> end do    ; end of station loop
> end do    ; end of time loop
>
>     npts=ntimes
>     fName =
> ip+"T2oben_1kmtemp_rainnonc_rainc_raintot_wspd_wdir_snow_snowh.txt"
>     data  = new( npts, "string")
>     ;print(" Time temp ")
>    do  ip = 0, 18
>     do it = 0,ntimes-1 ;ntimes-1
>     data (it)= (sprintf("%5.0f",it)    +" " \
>           +sprintf("%21.2f", temp(it)) +"  " \
>           +sprintf("%22.2f", rain_expp(it)) +"  " \
>           +sprintf("%23.2f", rain_conn(it)) +"  " \
>           +sprintf("%24.2f", rain_tot(it)) +"  " \
>           +sprintf("%17.2f", wspd(it)) +"  " \
>           +sprintf("%25.2f", wdir(it)) +"  " \
>           +sprintf("%18.2f", snoww_eq(it)) +"  " \
>           +sprintf("%19.2f", snowh_h(it)) +"  " )
> end do
>  end do                     ; end of time loop
>
>  asciiwrite (fName , data)
>
> 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/20160509/a1d9d782/attachment.html 


More information about the ncl-talk mailing list