[ncl-talk] Average Spatial Correlations
Dennis Shea
shea at ucar.edu
Mon Jun 13 09:34:21 MDT 2016
I must leave so not much time to respond.
Correlation (r) is a bounded quantity: -1 to +1.
[1]
I think the correct way to compute a mean correlation is via the Fischer
z-transform
rz = 0.5*log((1+r)/(1-r)) ; Fischer z-transform of correlations
[2]
Average the z-transformed quantities: rz_avg
[3]
Perform the inverse transform.
r_avg = tanh(rz_avg) ; same as =(exp(2*rz)-1)/(exp(2*rz)+1)
Good luck
On Mon, Jun 13, 2016 at 9:15 AM, Melissa Lazenby <M.Lazenby at sussex.ac.uk>
wrote:
> Dear NCL Users
>
> I am wanting to create 1 plot of averaged spatial correlations.
> I currently have 40 global spatial correlation plots but would like to
> average those all out to 1 global plot with the overall average spatial
> correlations at each gridpoint.
>
> I am not sure whether to use the average function or to sum them up first
> and divide by 40.
>
> Your help would be very much appreciated.
>
> Below is the code I am trying to create this average plot from.
>
> Kindest Regards
> Melissa
>
> ; ==============================================================
> ; eof_1.ncl
> ; ==============================================================
>
> 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
> ; ==============================================================
> ; User defined parameters that specify region of globe and
> ; ==============================================================
> latS = -40.
> latN = 0.
> lonL = 10.
> lonR = 60.
>
> yrStrt = 1975
> yrLast = 2004
>
> season = "DJF" ; choose Dec-Jan-Feb seasonal mean
>
> neof = 1 ; number of EOFs
> optEOF = True
> optEOF at jopt = 0 ; This is the default; most commonly used; no need to
> specify.
> optETS = False
>
> model =
> (/"ACCESS1-0","ACCESS1-3","bcc-csm1-1","bcc-csm1-1-m","BNU-ESM","CanESM2","CCSM4","CESM1-BGC","CESM1-CAM5","CESM1-FASTCHEM","CESM1-WACCM","CMCC-CESM","CMCC-CM","CMCC-CMS","CNRM-CM5","CSIRO-Mk3-6-0","EC-EARTH","FGOALS-g2","GFDL-CM3","GFDL-ESM2G","GFDL-ESM2M","GISS-E2-H","GISS-E2-H-CC","GISS-E2-R","GISS-E2-R-CC","HadGEM2-AO","HadGEM2-CC","HadGEM2-ES","inmcm4","IPSL-CM5A-LR","IPSL-CM5A-MR","IPSL-CM5B-LR","MIROC5","MIROC-ESM","MIROC-ESM-CHEM","MPI-ESM-LR","MPI-ESM-MR","MPI-ESM-P","MRI-CGCM3","NorESM1-M","NorESM1-ME"/)
>
>
> wks = gsn_open_wks("X11","ALL_40_MMM_Eof_cor_SST") ;
> specifies a plot
> gsn_define_colormap(wks,"ncl_default") ; choose color map
> plot = new (dimsizes(model),"graphic")
> do gg = 0,dimsizes(model)-1
>
> ; read in model data
>
>
> in=addfile("/mnt/geog/ml382/melphd/eof_sicz/eof90models/pr_"+model(gg)+"_eof_DJF_Anom.nc","r")
>
>
> lat = in->lat
> lon = in->lon
> time = in->time
> YYYY = cd_calendar(time,-1)/100 ; entire file
> iYYYY = ind(YYYY.ge.yrStrt .and. YYYY.le.yrLast)
>
> if((gg.eq.19) .or. (gg.eq.20)) then ;(gg.eq.17) .or.
> PR = in->pr(iYYYY,:,:)
> else
> PR = in->pr(iYYYY,0,:,:)
> end if
> ;
> nyrs = dimsizes(PR&time)
> ; =================================================================
> ; normalize data at each gridpoint by local standard deviation at each
> grid pt
> ; =================================================================
> PR = dim_standardize_n(PR,1,0)
>
> ; =================================================================
> ; Reorder (lat,lon,time) the *weighted* input data
> ; Access the area of interest via coordinate subscripting
> ; =================================================================
> x = PR({lat|latS:latN},{lon|lonL:lonR},time|:)
>
> eof = eofunc_Wrap(x, neof, optEOF)
> eof_ts = eofunc_ts_Wrap (x, eof, optEOF)
>
> printVarSummary( eof ) ; examine EOF variables
> printVarSummary( eof_ts )
>
> ; =================================================================
> ; Extract the YYYYMM from the time coordinate
> ; associated with eof_ts [same as SLP&time]
> ; =================================================================
>
> yyyymm = cd_calendar(eof_ts&time,-2)/100
>
> ; ==============================================================
> ; Open the file: Read only the user specified period
> ; ==============================================================
> s = addfile
> ("/mnt/nfs2/geog/ml382/melphd/global/sstmodel90/tos_Omon_"+model(gg)+"_historical_safrica_climDJF.nc",
> "r")
>
>
> TIME1 = s->time
> YYYYZ = cd_calendar(TIME1,-1)/100 ; entire file
> iYYYYZ = ind(YYYYZ.ge.yrStrt .and. YYYYZ.le.yrLast)
>
> sst = s->tos(iYYYYZ,:,:)
> printVarSummary(sst)
>
> sst1 =sst(:,:,:)
>
>
> printVarSummary(eof_ts) ; variable overview
> eof1 = eof_ts(0,:)
>
> sst1 = lonFlip(sst1)
> ; =================================================================
> ; Reorder (latitude,longitude,time) the *weighted* input data
> ; Access the area of interest via coordinate subscripting
> ; =================================================================
>
> q = eof_ts(0,:)
>
> y = sst1(lat|:,lon|:,time|:)
>
> y&lat at units = "degrees_north"
> y&lon at units = "degrees_east"
>
> x&lat at units = "degrees_north"
> x&lon at units = "degrees_east"
>
> printVarSummary(q)
> printVarSummary(y)
>
>
> ccr = escorc(q, y)
>
> printVarSummary(ccr)
>
> ccr!0 = "lat" ; name dimensions
> ccr!1 = "lon"
> ccr&lat = y&lat ; assign coordinate values and
> ccr&lon = y&lon ; units attributes
>
> printVarSummary(ccr)
>
>
> ;ensccrdimavg = dim_avg_Wrap( ccr(lat|:, lon|:)
> ensccr = new((/dimsizes(lat),dimsizes(lon)/), typeof(ccr))
> ;ensccr!0="lat"
> ;ensccr&lat=lat
> ;ensccr!1="lon"
> ;ensccr&lon=lon
> ;ensccr(gg) = ensccrsum
>
> ;ensccrsum(gg) = sum( ccr(lat|:, lon|: ))
> ;ensccr = ensccrsum / 40
>
> ensccr = dim_avg_n_Wrap(ccr,0)
>
>
> ensccr!0 = "lat" ; name dimensions
> ensccr!1 = "lon"
> ensccr&lat = lat ; assign coordinate values and
> ensccr&lon = lon ; units attributes
>
> printVarSummary(ensccr)
> ;============================================================
> ; PLOT
> ;============================================================
>
>
> ;*******************************************
> ; first plot
> ;*******************************************
>
> rescn = True
> rescn at cnFillOn = True
> rescn at cnLinesOn = False ; True is default
> rescn at gsnDraw = False ; don't draw
> rescn at gsnFrame = False ; don't advance frame
> rescn at cnLineLabelsOn = False ; True is default
> rescn at lbLabelBarOn = False ; turn off individual lb's
>
> rescn at cnLevelSelectionMode = "ManualLevels" ; set manual contour
> levels
> rescn at cnMinLevelValF = -0.5 ; set min contour level
> rescn at cnMaxLevelValF = 0.5 ; set max contour level
> rescn at cnLevelSpacingF = 0.1 ; set contour spacing
>
> ;rescn at lbOrientation = "Vertical" ; vertical label bar
> rescn at tiMainString = ""+model(gg)+""
> ;---This resource defaults to True in NCL V6.1.0
> rescn at lbLabelAutoStride = True ; optimal label stride
>
> rescn at gsnSpreadColors = True ; use full range of
> colors
> rescn at mpCenterLonF = 180. ; center plot at 180
> rescn at lbLabelFontHeightF = 0.015
>
> rescn at gsnAddCyclic = True
> ; reverse the first two colors
> setvalues wks
> "wkColorMap" : "ncl_default"
> "wkForegroundColor" : (/0.,0.,0./)
> "wkBackgroundColor" : (/1.,1.,1./)
> end setvalues
>
>
> end do
>
> plot = gsn_csm_contour_map(wks,ensccr,rescn)
>
>
> draw(plot)
> frame(wks)
>
> 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/20160613/498057b1/attachment.html
More information about the ncl-talk
mailing list