Re: EOFs for differently dimensioned arrays (fwd)

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Fri, 04 May 2007 12:21:14 -0600

Before I know how to answer, I have to know what you are doing :-)
I am probably missing something.

[1] Given:
       x = fi->$varNam$(time|st:en, {level|plev}, {lat|minlat:maxlat},
{lon|minlon:maxlon})
       Hence: x(time,lat,lon)

[2] Then:
        do nl=0,nlat-1 ; wgt gph by
sqrt(cos(latitude))
            x(nl,:,:) = x(nl,:,:)*cos(lat(nl)*rad)
        end do

       does not make sense. You are looping over the leftmost dimension
which is "time",
       yet you are treating it like latitude. I think you want the
following. Also,
       I have included the sqrt(cos(lat)) to conform to the comment

       do nl=0,nlat-1 ; wgt gph by sqrt(cos(latitude))
           ; x(nl,:,:) = x(nl,:,:)*sqrt( cos(lat(nl)*rad) )
             x(:,nl,:) = x(:,nl,:)*sqrt( cos(lat(nl)*rad) )
       end do

      Note: http://www.ncl.ucar.edu/Document/Functions/Built-in/eofunc.shtml
                Example 7 has a more efficient though less clear way of
applting the weighting.

[3] check_3d = eofunc(x, 5, False)

      The above may not be doing what you think....
      The classic application would be check_3d =
eofunc(x(lat|:,lon|:,time|:), 5, False)
      The covariance matrix would contain the covariance between grid
points over time.

     Your usage s equivalent to check_3d =
eofunc(x(time|:,lat|:,lon|:), 5, False)
      The covariance matrix would contain the covariance between each
(time,lat) over lon.

     I'll stop here because I am confused. You can email me offline if
you want.

Regards
D

   
      Here

>
> I am performing a PC analysis on the SH NCEP 500hPa geopotential height
> field. For various reasons, I have converted the data from a 3d (time,
> lat, lon) array to a 2d array (time, locn). When I perform the PC
> analysis on the 2d array, the percentage variance explained by the first
> EOF is far too high, approx. 98% (should be nearer 30%, perhaps a bit
> higher since I haven't removed the annual cycle). I cannot think of a
> fundamental reason for this, so presumably there is a problem with how
> the array transformation is performed, but I have no idea what I'm doing
> wrong. The code is below...
>
>
>
> ;=====================================================================;
> 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/csm/shea_util.ncl"
> ;======================================================================
>
> begin
>
> ;******************************************
> ;Define hyperslab
> ;*****************************************
>
> maxlat = 90.
> minlat = -90.
>
> minlon = 0.
> maxlon = 360.
>
> plev = 500.
>
> yrst = 1980
> yren = 2004
>
> diri = "~/datasets/atm/"
> fili = "ncep_geoz_1960_2004"
> fyrst = 1960
>
> varNam = "hgt"
>
> ;********************************
> ;read in data
> ;*****************************
>
> ;start/end time index
>
> st = (yrst - fyrst) * 12
> en = ((yren - fyrst) * 12) + 11
>
> fi = addfile(diri+fili+".nc", "r")
>
>
>
> x = fi->$varNam$(time|st:en, {level|plev}, {lat|minlat:maxlat},
> {lon|minlon:maxlon})
>
>
> dimz = dimsizes(x)
> ntim = dimz(0)
> nlat = dimz(1)
> nloc = nlat * dimz(2)
>
>
> ; COSINE WGT
>
> lat = x&lat
>
> rad = 3.14157/180.
>
> do nl=0,nlat-1 ; wgt gph by
> sqrt(cos(latitude))
> x(nl,:,:) = x(nl,:,:)*cos(lat(nl)*rad)
> end do
>
>
>
> ;create 2-d array: rows = time, column = locn.
>
> x2d = new((/ntim, nloc/), typeof(x))
>
> do tt = 0, ntim-1
>
> x2d(tt,:) = ndtooned(x(tt,:,:))
>
> end do
>
> ;alternatively, I tried the method here with the same result...
> ; x2d = onedtond(ndtooned(x), (/ntim, nloc/))
>
> ;eof
>
>
> ;2d EOF - results in pcvar = 98%
>
> check_2d = eofunc(x2d, 5, False)
> print("2d eof - "+check_2d_at_pcvar)
>
> ;EOF on original data, pcvar is reasonable
>
> check_3d = eofunc(x, 5, False)
> print("3d eof - "+check_3d_at_pcvar)
>
>
> ;Rearrange 2d array back into 3ds, pcvar is same as for original 3d
> array case
> x3d = new(dimz, typeof(x))
>
> do tt = 0, ntim-1
> x3d(tt,:,:) = onedtond(x2d(tt,:),(/dimz(1), dimz(2)/))
> end do
>
> check_ret = eofunc(x3d, 5, opt)
> print("regrid "+check_ret_at_pcvar)
>
> end
>
>
>
>
> ============================
> Room A113
> UCLA Department of Geography
> 1255 Bunche Hall
> Los Angeles CA 90095-1524
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk_at_ucar.edu
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

-- 
======================================================
Dennis J. Shea                  tel: 303-497-1361    |
P.O. Box 3000                   fax: 303-497-1333    |
Climate Analysis Section                             |
Climate & Global Dynamics Div.                       |
National Center for Atmospheric Research             |
Boulder, CO  80307                                   |
USA                        email: shea 'at' ucar.edu |
======================================================
_______________________________________________
ncl-talk mailing list
ncl-talk_at_ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri May 04 2007 - 12:21:14 MDT

This archive was generated by hypermail 2.2.0 : Tue May 08 2007 - 09:54:33 MDT