[ncl-talk] Probability Calculation
Ipshita Majhi
ipmajhi at alaska.edu
Thu Aug 27 12:07:27 MDT 2015
Dear NCL,
I have been working on creating a plot with correlation and probability
function. I am able to do it with no lag between two variables, but when I
introduce lag, the probability function is not working, I am not sure why
that is, it will be great if you could guide me on this:-
;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")
;==============================================
;Copying Attributes and Coordinates
copy_VarAtts(sst_1967_2012,sst_djf)
copy_VarCoords_1(sst_1967_2012,sst_djf)
;Calculating correlation
djf_reorder=sst_djf(lat|:,lon|:,time|:)
;Calculating Correlation Coefficient
;Months with Lag
corr_djf=esccr(JJAS_1967_2012,djf_reorder,1)
;========================================
;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/20150827/a706cd3b/attachment.html
More information about the ncl-talk
mailing list