[ncl-talk] Calculate percentile
Dave Allured - NOAA Affiliate
dave.allured at noaa.gov
Tue Mar 15 09:50:27 MDT 2022
Georgio, you said you want to compute the 70th percentile month by month.
I recommend a single outer loop over time, and do all the percentile
calculations independently inside the time loop. Sort the data for only
one month at a time, excluding missing values as needed; get the 70th
percent index within that month; and that gives you the 70th percentile for
the current month.
As Jonathan said, you need to use ndtooned() to reshape the input data to
1-D for input to NCL's sort routines. Do this independently for each
month, inside the time loop.
A long time ago, NCL's dim_pqsort function was faster than the qsort
function, due to internal coding. If your data set is large, you might
want to check on which function is currently faster.
On Tue, Mar 15, 2022 at 7:02 AM Buzan, Jonathan via ncl-talk <
ncl-talk at mailman.ucar.edu> wrote:
> Watchout, there’s a typo. The indexes need to match the lat and lon.
>
> mastervar(j,k) = work( min((/nwork-1, nwork*70/100-1 /)) ) ; 70%
> Change to:
> mastervar(k,j) = work( min((/nwork-1, nwork*70/100-1 /)) ) ; 70%
>
> On Mar 15, 2022, at 13:55, Buzan, Jonathan via ncl-talk <
> ncl-talk at mailman.ucar.edu> wrote:
>
> Hi Giorgio,
>
> The code for the stat_dispersion is
> in $NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl
> (Search for stat_dispersion for the explicit function)
>
> Here are some snippets to help. Ultimately, the stat_dispersion operates
> on 1d fields, and you have to put your data into that format. Highly
> recommend that you have more than 100 data points in your time dimension or
> the function will not work.
>
> mastervar = ; new array that is the stat dispersion by lat by lon
> do k = 0, dimsizes(lat)-1
> do j = 0, dimsizes(lon)-1
> z = data(:,k,j)
> work = ndtooned(z)
> nwork = dimsizes(work)
> qsort(work) ; reordered from lowest to highest
> mastervar(j,k) = work( min((/nwork-1, nwork*70/100-1 /)) ) ; 70%
> delete(z)
> delete(work)
> delete(nwork)
> end do
> end do
>
> Cheers,
> -Jonathan
>
> On Mar 15, 2022, at 13:25, Giorgio Graffino via ncl-talk <
> ncl-talk at mailman.ucar.edu> wrote:
>
> Dear NCLers,
>
> I have a time,lat,lon monthly-mean SST field and I want to compute the
> 70th percentile month by month. I'm trying to adapt Dennis's suggestion given
> here (https://www.ncl.ucar.edu/Support/talk_archives/2013/0954.html), but
> I'm ending up with a single 2-D percentile field, calculated across the
> entire time period. The problem is that, since it's historical SST, most
> of the spatial domain is above the threshold at the end of the time
> period, because of global warming.
>
> I'd like to have something similar to stat_dispersion (
> https://www.ncl.ucar.edu/Document/Functions/Contributed/stat_dispersion.shtml),
> but for a custom percentile and without detrending the data. Is there a
> way to do that? Here is a snippet of my code:
>
> temp_sort = dim_pqsort_n(temp,2,0) ; (time, lat, lon) -> (time, lat, lon)
>
> temp_num = dim_num_n(.not.ismissing(temp_sort),0) ; (time, lat, lon) ->
> (lat, lon)
>
> do n = 0, ntim-1
>
> temp_warmest30(n,:,:) =
> where(temp_sort(n,:,:).gt.floor(0.7*temp_num),temp(n,:,:),temp at _FillValue)
> ; (time, lat, lon) -> (time, lat, lon)
>
> end do
>
> Thanks,
>
> Giorgio
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at mailman.ucar.edu
> List instructions, subscriber options, unsubscribe:
> https://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20220315/bfd2b44d/attachment.html>
More information about the ncl-talk
mailing list