[ncl-talk] calculating and plotting change in standard deviation and change of both skewness and kutosis at 99% percentile

[Sri.durgesh Nandini-Weiss]

Hello everyone,

I'm attempting to calculating, and plot changes for extreme value from 
sea surface height from an ensemble outputs.

For e.g i am trying to analyse the change between 1986-2006 and 
2081-2100 at the 99% percentiles from the ensemble from the change in 
statistics (mean,  standard deviation and  change of both skewness and 
kutosis). I adapted my script from extval_6.ncl

Is the below the correct way?

;************************************************
; extvalv_6.ncl
;
; Concepts illustrated:
;   - Using extval_mlegev  and extval_gev
;   - Calculate basic statistics
;     - Basic statistics of the original sample
;     - Estimate GEV distribution parameters
;   - Creating a 'text' object to attach to a plot
;************************************************

    fn  = "zo_hist_1850-2005_ens_1-100.nc" ; define filename [time | 
1872] x [ens | 100] x [depth | 1] x [lat | 45] x [lon | 90]

    in  = addfile(fn,"r")                                 ; open netcdf file
    ts  = in->zo                                          ; read in variable
    printVarSummary(ts)

;************************************************
;  calculate deviation from zonal mean
;************************************************
    anom = dim_rmvmean_n_Wrap(ts,0)
    printVarSummary(anom)

;***************************************************************
;--- "Normal" (conventional) statistical estimates for full sample
;***************************************************************

    xStat4 = dim_stat4_n(ts, 0)  ; 1st 4 moments of original sample
    printVarSummary(xStat4 )

    xAvg   = xStat4(0)              ; original sample mean
    xStd   = sqrt(xStat4(1))        ;    "     sample std dev
    xSkew  = xStat4(2)              ; skewness; departure from symmetry
    xKurt  = xStat4(3)              ; kurtosis; relative to a normal 
distribution
    xLow   = min(ts)
    xMed   = dim_median_n(ts,0)  ; median of original sample
    xHi    = max(ts)


;***************************************************************
;--- PLOTS
;***************************************************************

    wks  = gsn_open_wks ("pdf","ext6")                ; send graphics to 
PNG file

    gsn_define_colormap(wks,"default")
    plot = new(2, "graphic")

;************************************************
; Panel
;************************************************

    resP                     = True                  ; modify the panel plot
    resP at gsnMaximize         = True                  ; ps, eps, pdf
    resP at gsnPanelMainString  = "River Flow Rate"     ; use this for NCL 
V6.4.0 and later
    resP at txFontHeightF       = 0.020
    gsn_panel(wks,plot,(/1,2/), resP)                ; now draw as one plot

;***************************************************************
;--- create histogram for the original sample
;***************************************************************

    resh = True
    resh at gsnDraw      = False
    resh at gsnFrame     = False
    resh2gsnHistogramNumberOfBins = 11

    resh at gsFillColor  = "green"
    resh at tiMainString = "SSH: N="+ntim
    plt_hist          = gsn_histogram(wks, ts ,resh)

;***************************************************************
;--- text object original sample statistics; place on histogram
;***************************************************************

    txres                       = True
    txres at txFont                = "helvetica-bold"
    txres at txFontHeightF         = 0.0150
    textSample  = (/"  Mean="+sprintf("%5.1f", xAvg) +"~C~"+  \
                    "   Std="+sprintf("%5.1f", xStd) +"~C~"+  \
                    "  Skew="+sprintf("%5.2f", xSkew) +"~C~"+ \
                    "  Kurt="+sprintf("%5.2f", xKurt) +"~C~"+ \
                    "  xLow="+sprintf("%5.1f", xLow) +"~C~"+  \
                    "  xMed="+sprintf("%5.1f", xMed) +"~C~"+  \
                    "   xHi="+sprintf("%5.1f", xHi ) /)
    txBoxSample  = gsn_create_text(wks,textSample, txres)

    amres  = True
    amres at amParallelPosF   =  0.30             ; move legend to the right
    amres at amOrthogonalPosF = -0.10             ; move the legend up
    annoSample = gsn_add_annotation(plt_hist, txBoxSample, amres)  ; 
Attach string to plot

    draw(plt_hist)
    frame(wks)



the change in statistics-- 
  Sri Nandini-Weiß
Research Scientist

Universität Hamburg
Center for Earth System Research and Sustainability (CEN)
Cluster of Excellence 'Climate, Climatic Change, and Society' (CLICCS)

Bundesstrasse 53, 20146 Hamburg
Tel: +49 (0) 40 42838 7472

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20191030/bad08bee/attachment.html>