[ncl-talk] Lanczos Filtering

[Dave Allured - NOAA Affiliate]

Guigma,

Here is another idea.  Instead of filtering precipitation, consider
filtering log precipitation, then get the antilog.  The Lanczos or other
similar filter should be okay.  Depending on the rest of your application,
you may want to apply a proportional correction so that the mean output is
the same as the mean input over time.

I have not tested this filter method.  However, it is guaranteed to not
produce negative output values.  It may also be more mathematically valid
for the one-sided and spiky behavior of precipitation.

However, be careful of changes in certain derived statistics.  I suspect
that number of days with non-zero precip will probably increase
dramatically, and this would be an artificial and meaningless side effect
of filter application.  Similarly I would distrust any derivations of
statistics related to extreme values.  This leads me back to something I
said before, search to see how others have dealt with filtering
precipitation.

--Dave


On Tue, Dec 5, 2017 at 1:41 PM, Dave Allured - NOAA Affiliate <
dave.allured at noaa.gov> wrote:

> Some negative weights are normal for Lanczos filters and similar digital
> filters.  Filters with sharp cutoffs (narrow transition bands) will produce
> overshoot or "ringing" as a normal side effect.  Please see this article
> about overshoot.
>
>     https://en.wikipedia.org/wiki/Overshoot_(signal)
>
> Precipitation is badly behaved for the purpose of filtering because (1) it
> has a one sided distribution around and including zero, and (2) it is
> naturally extremely spiky.  These two factors conspire to create lots of
> overshoot and small negative values when you apply ordinary weighted
> average filters.
>
> I do not have a simple remedy for your precip filtering application.  Try
> some simple searches such as "filter precipitation time series" to see how
> others may have approached this.  As Dennis said, good luck.   ;-)
>
> --Dave
>
>
> On Tue, Dec 5, 2017 at 1:19 PM, Dennis Shea <shea at ucar.edu> wrote:
>
>> Likely,  you have zero  or near-zero values in the original time series.
>>
>> This make filtering problematical. If you have any negatives weights you
>> may have filtered values less than 0,0.
>>
>> The reason I pointed you to Example 8 was so you could see the response
>> function of your parameters you set..
>>
>> Run the attached script. It produces 2 plots one for each 'sigma' factor.
>>
>> You may have to resort to
>>
>>    prc_filt = where(prc_filt.lt.0.0, 0.0, prc_filt)
>>
>> Good luck
>>
>>
>> On Tue, Dec 5, 2017 at 12:37 PM, Kiswendsida Hyacinthe GUIGMA <
>> karongseba at gmail.com> wrote:
>>
>>> Thank you dear Dennis for your help. I am sorry if my questions look too
>>> basic. I already read the filtering examples.
>>> The fact is that using their method (you can see my code in my previous
>>> post), from daily rainfall data I got negative values after filtering, what
>>> is abnormal.
>>> -----------------------------
>>> GUIGMA
>>>
>>>
>>> On Tuesday, 5 December 2017, 19:03, Dennis Shea <shea at ucar.edu> wrote:
>>>
>>> Please look at the filtering examples:
>>>    https://www.ncl.ucar.edu/Applications/filter.shtml
>>>
>>>    In particular, run Example 8
>>>
>>> [1] yes, of course ... more weights can create a sharper cutoff but
>>> there is a larger loass od data at the beginning and end.
>>>
>>> [2] re: practical criterion .... That is up to you as a scientist.
>>>
>>> --
>>> Note ncl-talk is not the place to learn about filtering.
>>>
>>> Good luck
>>>
>>>
>>> On Tue, Dec 5, 2017 at 9:14 AM, Kiswendsida Hyacinthe GUIGMA <
>>> karongseba at gmail.com> wrote:
>>>
>>> Dear NCL users,
>>> I am very new on NCL and I am trying to filter my rainfall data using
>>> Lanczos Filter (bandpass 10-90 day filtering). My question concerns the
>>> processing of the weights. I would like to know
>>> 1-is there any relationship between the number of weights to choose and
>>> the cutoff frequencies?
>>> 2-is there any practical criterion to choose the number of weights? I
>>> know that this number influences the loss of data.
>>> Below is my code:
>>>
>>> ------------------------------ ------------------------------
>>> ------------------------------ ------------------------------
>>> ------------------------------ ----------------
>>> ;import the file
>>> gpcp=addfile("gpcp_Sahel.nc"," r")
>>> rr=gpcp->precip
>>> ;filter parameters
>>> ihp=2 ;bandpass
>>> sigma=1.0 ; Lanczos sigma
>>> nWgt    = 81  ; number of weights -> means 40 data lost at the beginning
>>> and 40 data lost at the end.
>>> fca     = 1./90.            ; start freq ;highpass= 90 days
>>> fcb     = 1./10.                        ; last  freq; lowpass= 10 days
>>> wgt     = filwgts_lanczos (nWgt, ihp, fca, fcb, sigma )
>>> printVarSummary(wgt)
>>> xBPF    = wgt_runave_n_Wrap ( rr, wgt, 0, 0)      ; 10-90 day bp
>>> performed across time dimension(0) of rr (meaning the time dimension)
>>>  ----------------------------- ------------------------------
>>> ------------------------------ ------------------------------
>>> ------------------------------ ----------------
>>>
>>> Cheers,
>>>
>>> GUIGMA
>>>
>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171205/199f8f8b/attachment.html>