Issue with sortl: [ferret_users] Percentile along time axis for gridded data

I am using the method provided by Russ a couple of years ago to compute percentile along time axis for gridded data of sea ice concentrations.

The resultant Q-tiles are a bit of too scattered in the seasonal varying ice zone. I have worked with Russ directly during the shutdown – Thank you, Russ, for your helpful suggestions.

looking into it in more details, I am a bit of confused by the behavior of function sortl which does not seem to sort correctly in ldim at each I,j grid cell in the sense that the largest value should be sorted as the last. It can be demonstrated by the following ferret output. (I was using a much earlier version but have reinstalled the latest per Russ’s suggestion. The results are the same.)

Is this correct?

Thanks,

--- Peng

$ ferret

NOAA/PMEL TMAP

FERRET v7.44 (optimized)

Linux 3.10.0-957.1.3.el7.x86_64 64-bit - 12/07/18

4-Feb-19 14:01

yes? use coads_climatology

yes? let numpoints=`sst,return=lend`

!-> DEFINE VARIABLE numpoints=12

yes? let numlats=`sst,return=jend]`

!-> DEFINE VARIABLE numlats=90

yes? let/title="sorted indices" sorted_by_time=sortl(sst[l=1:`numpoints`])

!-> DEFINE VARIABLE/title="sorted indices" sorted_by_time=sortl(sst[l=1:12])

yes? list sst[i=90,j=45],sorted_by_time[i=90,j=45]

WARNING: Listed variables have ambiguous coordinates on axes: T

DATA SET: /snfs1/users/gpeng/myFerret/fer_dsets/FerretDatasets-7.4/data/coads_climatology.cdf

LONGITUDE: 161W

LATITUDE: 1S

Column 1: SST[T=01-JAN 00:45:31-DEC 06:34] is SEA SURFACE TEMPERATURE (Deg C)

Column 2: SORTED_BY_TIME[T=0.5:12.5] is sorted indices

SST SORTED_BY_TIME

L / 1: 26.82 1.00

L / 2: 26.99 2.00

L / 3: 27.47 12.00

L / 4: 27.83 11.00

L / 5: 27.96 9.00

L / 6: 28.22 8.00

L / 7: 27.96 3.00

L / 8: 27.42 10.00

L / 9: 27.33 4.00

L / 10: 27.55 5.00

L / 11: 27.22 7.00

L / 12: 27.02 6.00

On Tue, Mar 29, 2016 at 10:16 PM Russ Fiedler <russell.fiedler@xxxxxxxx> wrote:

Hi Peng,

I did this in the dim distant past by saving the sorted data to a file and then picked the appropriate values out.

use coads_climatology
let numpoints=`sst,return=lend`
let numlats=`sst,return=jend]`
let/title="sorted indices" sorted_by_time=sortl(sst[l=1:`numpoints`])

! I found that I ran out of memory if I tried to sort too much for some GB sized datasets. It depends on the size of your dataset.
! You might be able to make larger windows or even do it in one go.

save/j=1/jlimits=1:`numlats`/clob/file=sorted.nc sorted_by_time
repeat/j=2:`numlats` save/app/file=sorted.nc sorted_by_time

can var sorted_by_time
use sorted

! Make a variable with time indices

let index_2d=0*x[g=sorted_by_time] + 0*y[g=sorted_by_time] + l[g=sorted_by_time]

! Create an integrating kernel for the percentile wanted 1 for the limit we want and missing elsewhere.

!25 percentile say. Make sure at least 1 valid point is used

let k25 = if l[g=index_2d] eq max(int(0.25*ignore0(sorted_by_time[l=1:`numpoints`@ngd])),1) then 1

! Show the indices used in the sorted data. Note that missing values are accounted for.

shade k25[l=@loc:1]

! Multiply by SST and sum up

let sst_k25 = k25[gl=sst[d=coads_climatology]@asn]*sst[d=coads_climatology]

let/title="sst at 25 percentile" sst25=sst_k25[l=1:`numpoints`@sum]

shade sst25

Cheers,
Russ

On 30/03/16 05:32, Ge Peng - NOAA Affiliate wrote:

Found this message showing how to find percentiles for 2-dimensional gridded data:

http://ferret.pmel.noaa.gov/Ferret/Mail_Archives/fu_2004/msg00632.html

However, I would like to compute quantiles/percentiles along my time axis at each grid cell of the 2-dimensional gridded data. I.e., for each x and y, the percentiles are done along the time axis with valid data points. (We can ignore the z dimension for now.)

I could take the time series at each grid cell using nested repeat loop in x and y dimensions, following the above example to sort the data and place the ordered data onto the tiled axis.

It does not sound very efficient to me. Has anyone done something similar in ferret? Is there a better way to do this, perhaps without nested repeat loops in both x and y directions?

Appreciate any help.

--- Peng

--

Ge Peng, Ph.D
Research Scholar

Cooperative Institute for Climate and Satellites, NC (CICS-NC)

North Carolina State University (NCSU) and

NOAA’s National Centers for Environmental Information (NCEI)

Center for Weather and Climate (CWC)

151 Patton Ave, Asheville, NC 28801
ge.peng@xxxxxxxx
o: +1 828 257 3009
f: +1 828 257 3002

Following CICS-NC on Facebook

Ge Peng, PhD
Research Scholar
Cooperative Institute for Climate and Satellites - NC (CICS-NC)/NCSU at

NOAA’s National Centers for Environmental Information (NCEI)

Center for Weather and Climate (CWC)

151 Patton Ave, Asheville, NC 28801
+1 828 257 3009; ge.peng@xxxxxxxx

ORCID: http://orcid.org/0000-0002-1986-9115

Following CICS-NC on Facebook