[ncl-talk] extreme value analysis

[Sri.durgesh Nandini-Weiss]

Hello everyone


I tried calculating and plotting extreme value analysis on my SSH data, 
following ncl example extval_6.ncl 
<https://www.ncl.ucar.edu/Applications/Scripts/extval_6.ncl>:

My plot was generated with no error, but it doe not look correct.

See attached, as well as my script. Can someone advice me? I am happy to 
send my data as well.

Sri

;************************************************
; extvalv_6.ncl

    f     = addfile ("anom_hist_zo_slr.nc", "r")
    hist_anom    = f->hist_anom
    printVarSummary(hist_anom)


    dimx = dimsizes(hist_anom)
    ntim = dimx(0)           ; 240
    nens = dimx(1)           ; 100
    nlat = dimx(2)           ; 45
    nlon = dimx(3)           ; 90
    lat1=nlat
    lon1=nlon
    nmos = 12
    nyrs   = ntim/nmos       ; 20
    printVarSummary(nyrs)

;==================================================================
; print ("Data read in - start calculations")
; Region domain lat and lon
;==================================================================

   hist_anom1=hist_anom(:,:,{48},{-125})               ;(ntim, mlon)
   printVarSummary(hist_anom1)
;=================================================================
; extract data between 95-100percentiles
;=================================================================

   anom1_1d = ndtooned(hist_anom1)
   printVarSummary(anom1_1d)

   anom1_nomiss=anom1_1d(ind(anom1_1d.ne.-9999))
   anom1_permsort = dim_pqsort(anom1_nomiss,2)
   x95=anom1_nomiss(floattoint(0.95*dimsizes(anom1_nomiss)))
   i95  =  ind(anom1_nomiss.ge.x95)
   xHigh = anom1_nomiss(i95)                 ; values
   printVarSummary(xHigh)

;***************************************************************
;---GEV parameter estimation
;***************************************************************

    vals     = extval_mlegev (xHigh, 0, False)    ; dims=0
    print(vals)

    center   = vals(0)                            ; extract for clarity
    scale    = vals(1)
    shape    = vals(2)
    ;ptype=0
    pdfcdf   = extval_gev(xHigh, shape, scale, center,0)
    ;pdfcdf = extval_pareto(xHigh, shape, scale, center, ptype, 0)   ; 
[/ PDF, CDF /]
    pdf      = pdfcdf[0]                          ; extract from list 
for convenience
    cdf      = pdfcdf[1]
    delete(pdfcdf)

    printVarSummary(pdf)
    printMinMax(pdf,1)
    print("---")
    printVarSummary(cdf)
    printMinMax(cdf,1)
    print("---")

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

    plotfile="Hist_SSH"
   wks_type = "pdf"
   wks = gsn_open_wks(wks_type,plotfile)

    gsn_define_colormap(wks,"default")
    plot = new(2, "graphic")
    res     = True
    res at gsnDraw  = False
    res at gsnFrame = False
    res at tiYAxisString = ""

    it      = dim_pqsort_n(xHigh, 1, 0)        ; for plot reasons make  
'x' is asebding

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

;************************************************
;--- Plot PDF and CDF ; create panel
;************************************************

    plot(0) = gsn_csm_xy (wks,xHigh(it),pdf(it),res) ; create plot

    res at trYMinF  = 0.0
    res at trYMaxF  = 1.0

    plot(1) = gsn_csm_xy (wks,xHigh(it),cdf(it),res) ; create plot

    gsn_panel(wks,plot,(/1,2/), resP)             ; now draw as one plot


-- 
Dr. 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/20191212/4f635ac8/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Hist_SSH.pdf
Type: application/pdf
Size: 39382 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20191212/4f635ac8/attachment.pdf>