Re: Calculate time-averaged wind speed from 2 nc files

You can treat the 2 files as a single file using the "addfiles" function:

ncl 0> files = systemfunc("ls C4IRCA3_A2*")
ncl 1> f1 = addfiles(files,"r")
ncl 2> print(f1)

Variable: f1
Type: list
Total Size: 4 bytes
            1 values
Number of Dimensions: 1
Dimensions and sizes: [1]
Coordinates:
Type: list <concat | fifo>
Total items: 2

List Item 0: NclFileVarClass
Variable: unnamed
Type: file
File path: C4IRCA3_A2_ECHAM5_DM_25km_1981-1990_wss.nc
Number of global attributes: 8
Number of dimensions: 6
Number of variables: 9

List Item 1: NclFileVarClass
Variable: unnamed
Type: file
File path: C4IRCA3_A2_ECHAM5_DM_25km_1991-2000_wss.nc
Number of global attributes: 8
Number of dimensions: 6
Number of variables: 9

ncl 3> wss = f1[:]->wss
ncl 4> printVarSummary(wss)

Variable: wss
Type: float
Total Size: 1054842000 bytes
            263710500 values
Number of Dimensions: 4
Dimensions and sizes: [time | 7305] x [height | 1] x [rlat | 190] x [rlon | 190]
Coordinates:
            time: [11323.5..18627.5]
            height: [10..10]
            rlat: [-20.68..20.9]
            rlon: [-26.12..15.46]
Number Of Attributes: 7
  standard_name : wind_speed
  long_name : 10-meter wind speed
  units : m s-1
  cell_methods : time: mean
  coordinates : lon lat
  grid_mapping : rotated_pole
  _FillValue : 1e+30

ncl 5> printVarSummary(f1[0]->lat)

Variable: lat (file variable)
Type: float
Total Size: 144400 bytes
            36100 values
Number of Dimensions: 2
Dimensions and sizes: [rlat | 190] x [rlon | 190]
Coordinates:
            rlat: [-20.68..20.9]
            rlon: [-26.12..15.46]
Number Of Attributes: 3
  standard_name : latitude
  long_name : latitude
  units : degrees_north

Note that you use all the files ( f[:]->wss ) to get the data variable and the time variable, but just a single file to get the latitude ( f[0] ).
See the documentation on addfiles.
 -dave

On Apr 2, 2013, at 6:00 AM, Ioannis Koletsis <koletsis@noa.gr> wrote:

> 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
>
>
>
> <average3.ncl>_______________________________________________
> 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 Tue Apr 2 14:02:09 2013