[ncl-talk] Significance calculation

Ipshita Majhi ipmajhi at alaska.edu
Sat Aug 29 23:06:06 MDT 2015 

Dear NCL,

I would like to calculate the significance. I am using the probability
function, i am not sure where things are not working out. There are no
error message, the code runs but no stippling show for significance.
Here is the code

;This is for correlating JJAS monthly precp with JJAS temperature with lag

;******************************************************
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"
;*******************************************************

;=============================================
; Reading data
;=============================================

a=addfile("~/Documents/SST/Monthly/sst.mnmean.v4.nc","r")

;==============================================

JJAS_1871_2012
=asciiread("~/Documents/SST/Monthly/jjas_1871_2012.txt",(/142,1/),"float")
;Extracting JJAS for 1967 to 2012

JJAS_1967_2012=JJAS_1871_2012(96:141,0)
print(dimsizes(JJAS_1967_2012))

e=transpose(transpose(JJAS_1967_2012))
;print(dimsizes(e))

;=============================================

;For some months around september it is better to
;=============================================
 ;The dataset goes from 1854-2015
  sst=a->sst
; copy meta data and coordinate variables using contributed functions
 copy_VarAtts(a,sst)
 copy_VarCoords_1(a,sst)
;============================================

;Extracting 1967-2012 for SST
sst_1967_2012=sst(1367:1918,:,:)


copy_VarAtts(sst,sst_1967_2012)
copy_VarCoords_1(a,sst_1967_2012)

;============================================
;Seasonal three month average


 sst_djf=month_to_season(sst_1967_2012, "DJF")
 sst_jfm=month_to_season(sst_1967_2012, "JFM")
 sst_fma=month_to_season(sst_1967_2012, "FMA")
 sst_mam=month_to_season(sst_1967_2012, "MAM")
 sst_amj=month_to_season(sst_1967_2012, "AMJ")
 sst_mjj=month_to_season(sst_1967_2012, "MJJ")
 sst_jja=month_to_season(sst_1967_2012, "JJA")
 sst_jas=month_to_season(sst_1967_2012, "JAS")
 sst_aso=month_to_season(sst_1967_2012, "ASO")
 sst_son=month_to_season(sst_1967_2012, "SON")
 sst_ond=month_to_season(sst_1967_2012, "OND")
 sst_ndj=month_to_season(sst_1967_2012, "NDJ")

;==============================================
;Copying Attributes and Coordinates

copy_VarAtts(sst_1967_2012,sst_djf)
copy_VarCoords_1(sst_1967_2012,sst_djf)

copy_VarAtts(sst_1967_2012,sst_jfm)
copy_VarCoords_1(sst_1967_2012,sst_jfm)

copy_VarAtts(sst_1967_2012,sst_fma)
copy_VarCoords_1(sst_1967_2012,sst_fma)

copy_VarAtts(sst_1967_2012,sst_mam)
copy_VarCoords_1(sst_1967_2012,sst_mam)

copy_VarAtts(sst_1967_2012,sst_amj)
copy_VarCoords_1(sst_1967_2012,sst_amj)

copy_VarAtts(sst_1967_2012,sst_mjj)
copy_VarCoords_1(sst_1967_2012,sst_mjj)

copy_VarAtts(sst_1967_2012,sst_jja)
copy_VarCoords_1(sst_1967_2012,sst_jja)

copy_VarAtts(sst_1967_2012,sst_jas)
copy_VarCoords_1(sst_1967_2012,sst_jas)

copy_VarAtts(sst_1967_2012,sst_aso)
copy_VarCoords_1(sst_1967_2012,sst_aso)

copy_VarAtts(sst_1967_2012,sst_son)
copy_VarCoords_1(sst_1967_2012,sst_son)

copy_VarAtts(sst_1967_2012,sst_ond)
copy_VarCoords_1(sst_1967_2012,sst_ond)

copy_VarAtts(sst_1967_2012,sst_ndj)
copy_VarCoords_1(sst_1967_2012,sst_ndj)

;=============================================

;Calculating correlation



djf_reorder=sst_djf(lat|:,lon|:,time|:)
jfm_reorder=sst_jfm(lat|:,lon|:,time|:)
fma_reorder=sst_fma(lat|:,lon|:,time|:)
mam_reorder=sst_mam(lat|:,lon|:,time|:)
amj_reorder=sst_amj(lat|:,lon|:,time|:)
mjj_reorder=sst_mjj(lat|:,lon|:,time|:)
jja_reorder=sst_jja(lat|:,lon|:,time|:)
jas_reorder=sst_jas(lat|:,lon|:,time|:)
aso_reorder=sst_aso(lat|:,lon|:,time|:)
son_reorder=sst_son(lat|:,lon|:,time|:)
ond_reorder=sst_ond(lat|:,lon|:,time|:)
ndj_reorder=sst_ndj(lat|:,lon|:,time|:)

;===================================================

copy_VarAtts(sst_djf,djf_reorder)
copy_VarCoords_1(sst_djf,djf_reorder)

copy_VarAtts(sst_jfm,jfm_reorder)
copy_VarCoords_1(sst_jfm,jfm_reorder)

copy_VarAtts(sst_fma,fma_reorder)
copy_VarCoords_1(sst_fma,fma_reorder)

copy_VarAtts(sst_mam,mam_reorder)
copy_VarCoords_1(sst_mam,mam_reorder)

copy_VarAtts(sst_amj,amj_reorder)
copy_VarCoords_1(sst_amj,amj_reorder)

copy_VarAtts(sst_mjj,mjj_reorder)
copy_VarCoords_1(sst_mjj,mjj_reorder)

copy_VarAtts(sst_jja,jja_reorder)
copy_VarCoords_1(sst_jja,jja_reorder)

copy_VarAtts(sst_jas,jas_reorder)
copy_VarCoords_1(sst_jas,jas_reorder)

copy_VarAtts(sst_aso,aso_reorder)
copy_VarCoords_1(sst_aso,aso_reorder)

copy_VarAtts(sst_son,son_reorder)
copy_VarCoords_1(sst_son,son_reorder)

copy_VarAtts(sst_ond,ond_reorder)
copy_VarCoords_1(sst_ond,ond_reorder)

copy_VarAtts(sst_ndj,ndj_reorder)
copy_VarCoords_1(sst_ndj,ndj_reorder)

;=======================================
;Calculating Correlation Coefficient

;Months with Lag
corr_djf=esccr(JJAS_1967_2012,djf_reorder,1)
corr_ndj=escorc(ndj_reorder(:,:,0:44),JJAS_1967_2012(1:45))
corr_ond=escorc(ond_reorder(:,:,0:44),JJAS_1967_2012(1:45))
corr_son=escorc(son_reorder(:,:,0:44),JJAS_1967_2012(1:45))
corr_aso=escorc(aso_reorder(:,:,0:44),JJAS_1967_2012(1:45))

copy_VarAtts(djf_reorder,corr_djf)
copy_VarCoords_1(djf_reorder,corr_djf)

copy_VarAtts(ndj_reorder,corr_ndj)
copy_VarCoords_1(ndj_reorder,corr_ndj)

copy_VarAtts(ond_reorder,corr_ond)
copy_VarCoords_1(ond_reorder,corr_ond)

copy_VarAtts(son_reorder,corr_son)
copy_VarCoords_1(son_reorder,corr_son)

copy_VarAtts(aso_reorder,corr_aso)
copy_VarCoords_1(aso_reorder,corr_aso)

;========================================
;Probability
x=corr_djf(:,:,1)
print(dimsizes(x))
;==============================
; Now plotting the correlation

;===========================================

; Now plotting the correlation


    wks  = gsn_open_wks("pdf","DJF_sst_jjas_cor_prob_sig_lag")
  gsn_define_colormap(wks,"BlWhRe")              ; choose colormap

  res                      = True                ; make plot mods
  res at gsnDraw = False
  res at gsnFrame = False
  res at cnFillOn             = True                ; turn on color
  res at gsnSpreadColors      = True                ; use full colormap
  res at lbLabelAutoStride    = True         ; automatic lb label stride

  res at cnLinesOn            = False               ; turn off contour lines
  res at cnLevelSelectionMode = "ManualLevels" ; manually set cnlev
  res at cnMinLevelValF       = -1.                 ; min level
  res at cnMaxLevelValF       =  1.                 ; max level
  res at cnLevelSpacingF      = .1                  ; contour level spacing

  lag                      = 0
  res at tiMainString         = "DJF SST with JJAS rainfall"
  plotA = gsn_csm_contour_map_ce(wks,corr_djf,res)

 ; Significance plot on top of correlation plot
 ;========================= PLOT 2 =========
  res2 = True                            ; res2 probability plots

  res2 at gsnDraw             = False       ; Do not draw plot
  res2 at gsnFrame            = False       ; Do not advance frome
  res2 at cnFillOn = True                   ; turn on color/pattern fill
  res2 at cnMonoFillPattern = False   ; allow different patterns
  res2 at cnMonoFillColor = True       ; only use one color (black)

  res2 at cnLevelSelectionMode = "ExplicitLevels" ; set explicit cnlev
  res2 at cnLevels   = (/.90/)    ; only have 1 contour level
  res2 at cnFillPatterns = (/-1,17/) ; don't fill <0.95, stipple >=0.95

  res2 at gsnAddCyclic = True   ; may or may not be needed

  res2 at cnInfoLabelOn       = False       ; turn off info label

  res2 at cnLinesOn           = False       ; do not draw contour lines
  res2 at cnLineLabelsOn      = False    ; do not draw contour labels

  res2 at cnFillScaleF        = 0.5         ; add extra density

  printMinMax(prob_djf,0)   ; check to make sure values are
                                             ; between 0 and 1.

  plot2   = gsn_csm_contour(wks,prob_djf, res2)
  overlay (plotA, plot2)

  draw (plotA)
  frame(wks)

Best Regards
Ipshita
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150830/28645296/attachment.html

More information about the ncl-talk mailing list