Calculate time-averaged wind speed from 2 nc files

From: Ioannis Koletsis <koletsis_at_nyahnyahspammersnyahnyah>
Date: Tue Apr 02 2013 - 06:00:32 MDT

Dear NCL users,

I would like to calculate the average wind speed value of each grid point
for the period 1981-2000 from 2 .nc files, which have the same variables but
for different time periods.

(The .nc files are C4IRCA3_A2_ECHAM5_DM_25km_1981-1990_wss.nc and
C4IRCA3_A2_ECHAM5_DM_25km_1991-2000_wss.nc, available here:
http://ensemblesrt3.dmi.dk/data/A2/C4I/DM/)

 

Each of .nc files comprise of 10 year wind data (e.g. 1981-1990, 1991-2000)
and the attached script works fine for each file separately.

 

However, I am interested to calculate the average wind speed for the entire
period.

 

I tried to use a do loop before line 16 (do ifil=0, 1) in order to read the
files, but unfortunately does not work.

 

I think that merging the two files into one would be a solution, but I don't
know if there is an easy way in ncl..

 

Searching the ncl users list, I could not find any solution. Also, I would
like to get away from using NCO operators..

Any help would be appreciated..

 

John

 

 

 

SCRIPT

;*************************************************

; Define time-averaged wind field

;*************************************************

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"

 

begin

 

  DataDir = "/data3/koletsis/"

  Files = systemfunc (" ls -1 " + DataDir + \

         "C4IRCA3_A2_ECHAM5_DM_25km_* ")

  numFiles = dimsizes(Files)

  print(Files)

  print(" ")

 

  ifil = 0

  f1 = addfile(Files(ifil)+".nc","r")

  

 

  time = f1->time

  ws = f1->wss

  lat2d = f1->lat

  lon2d = f1->lon

 

 

  month_abbr = (/"","Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep", \

                    "Oct","Nov","Dec"/)

 

  time2 = ut_calendar(time, 0)

 

  year = tointeger(time2(:,0)) ; Convert to integer for

  month = tointeger(time2(:,1)) ; use sprinti

  day = tointeger(time2(:,2))

 

  date_str = sprinti("%0.2i ", day) + month_abbr(month) + " " + \

           sprinti("%0.4i ", year)

 

  beg_year = 1981

  end_year = 2000

  beg_month = 01

  end_month = 12

  beg_day = 01

  end_day = 31

 

  itime = ind(year.ge.beg_year.and.year.le.end_year.and.\

          month.ge.beg_month.and.month.le.end_month.and.\

          day.ge.beg_day.and.day.le.end_day)

 

; wsNew = ws(height|:,rlat|:,rlon|:,time|ibeg:iend)

; wsAvg = dim_avg_n_Wrap(ws,0)

  wsAvg = dim_avg_n_Wrap(ws(itime,:,:,:),0)

; print(wsNew)

; print(ibeg)

; print(iend)

   print(wsAvg)

; printVarSummary(wsAvg)

; asciiwrite("aver.txt", wsAvg)

 

end

 

 

 

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk

Received on Tue Apr 2 05:59:44 2013

This archive was generated by hypermail 2.1.8 : Tue Apr 02 2013 - 21:23:48 MDT