Re: Writing time series ascii file for each grid cell - Is there a faster way?

From: Dave Allured <dave.allured_at_nyahnyahspammersnyahnyah>
Date: Thu Jan 19 2012 - 16:26:01 MST

For starters, try switching from asciiwrite to write_matrix. See
example 6 on the write_matrix doc page. I believe that multiple "dat"
assignments and sprintf calls are extremely consumptive of CPU time
and memory in NCL, for large arrays; and this is within the innermost
loop where it really matters.

For write_matrix, you will need an intermediate array ("p" in example
6) that contains all the numeric data for the rows and columns to be
output. You will need to convert integers to floats for some of the
input data; date/time values in particular.

You can save more time and memory if you make this a 3-D array rather
than 2-D, and load it directly from the file input statements which
are in your outer loop, omitting the other input arrays pr, sph, srad,
etc. Also avoid deleting the intermediate array if you don't need to,
because it does not change shape between iterations.

--Dave

On Thu, Jan 19, 2012 at 12:16 PM, Bridget Thrasher
<bthrasher@climatecentral.org> wrote:
> I need to take many years of daily data for 6 variables and output their
> time series together to individual ASCII files for each of 273252 grid
> cells. I have code that works (see below), but when testing with a 32-year
> series, it's looking to take 3+ weeks to finish. I will soon need to do this
> many times over a 150-year series. Yikes! Am I going about this the wrong
> way?
>
> -Bridget
>
> begin
>
> ls_pr = systemfunc("ls 6km/pr/*.nc")
> ls_sph = systemfunc("ls 6km/sph/*.nc")
> ls_srad = systemfunc("ls 6km/srad/*.nc")
> ls_tmax = systemfunc("ls 6km/tasmax/*.nc")
> ls_tmin = systemfunc("ls 6km/tasmin/*.nc")
> ls_wind = systemfunc("ls 6km/wind/*.nc")
>
> f_pr = addfiles(ls_pr,"r")
> f_sph = addfiles(ls_sph,"r")
> f_srad = addfiles(ls_srad,"r")
> f_tmax = addfiles(ls_tmax,"r")
> f_tmin = addfiles(ls_tmin,"r")
> f_wind = addfiles(ls_wind,"r")
>
> lats = f_pr[0]->lat
> lons = f_pr[0]->lon
> times = f_pr[:]->time
> t = ut_calendar(times,-5)
> ts = dimsizes(times)
>
> do i = 0,dimsizes(lats)-1
> pr = f_pr[:]->pr(1:,i,:)
> sph = f_sph[:]->sph(1:,i,:)
> srad = f_srad[:]->srad(1:,i,:)
> tmax = f_tmax[:]->tasmax(1:,i,:)
> tmin = f_tmin[:]->tasmin(1:,i,:)
> wind = f_wind[:]->wind(1:,i,:)
> do j = 0,dimsizes(lons)-1
> if
> (ismissing(pr(0,j)).or.ismissing(sph(0,j)).or.ismissing(srad(0,j)).or.ismissing(tasmax(0,j)).or.ismissing(tasmin(0,j)).or.ismissing(wind(0,j)))
> then
> print("No data at "+lats(i)+"_"+(lons(j)-360.))
> else
> print("Making "+"data_"+lats(i)+"_"+(lons(j)-360.))
> fname = "ascii/data_"+lats(i)+"_"+(lons(j)-360.)
> dat = new((/ts-1/),"string","")
>
> dat = sprinti("%0.4i",t(1:,0))
> dat = dat + sprinti("%7.2i",t(1:,1))
> dat = dat + sprinti("%7.2i",t(1:,2))
> dat = dat + sprintf("%10.5f",pr(:,j))
> dat = dat + sprintf("%10.5f",tmax(:,j))
> dat = dat + sprintf("%10.5f",tmin(:,j))
> dat = dat + sprintf("%10.5f",wind(:,j))
> dat = dat + sprintf("%10.5f",sph(:,j))
> dat = dat + sprintf("%10.5f",srad(:,j))
>
> asciiwrite(fname,dat)
> delete(dat)
> end if
> end do
> delete(pr)
> delete(sph)
> delete(srad)
> delete(tmax)
> delete(tmin)
> delete(wind)
>
> end do
>
> end
>
>
> --
> Bridget Thrasher, PhD
> Independent Contractor, Research Scientist
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Jan 19 16:26:12 2012

This archive was generated by hypermail 2.1.8 : Mon Jan 23 2012 - 12:45:14 MST