Skip to content

Comments

New feature: Data detrending#435

Open
billetd wants to merge 7 commits intodevelopfrom
enh/detrend
Open

New feature: Data detrending#435
billetd wants to merge 7 commits intodevelopfrom
enh/detrend

Conversation

@billetd
Copy link
Contributor

@billetd billetd commented Feb 7, 2026

Scope

This PR adds functionality to detrend fitacf-level data with a running mean low-pass filter. This is useful for the detection of periodic fluctuations (e.g., ULF waves, TIDs), or for simply removing a background feature.

Currently, the detrending functionality works for velocities, SNR's, or both. I didn't see why it would be useful for spectral widths and elevation angles, but it would be trivial to add them if a reason emerged.

The code will work with any experiment where the beams are periodically repeating (e.g. normalscan, camping beam, widebeam, Themisscan). Care must be taken if used with irregular beam intervals.

It would also be pretty easy to slot in additional low-pass filter types in the future - just add another filter function to the class.

Calling sequence

dmap_detrended = Detrend.detrend_fitacf(fitacf_data, parameter, window_length, detrend_type)
"""
Detrend a series of input fitacf records using a low-pass moving average filter.
Useful for the study of period ULF waves or generally removing background flows.

Parameters
-----------
fitacf_data: List[dict]
    List of dictionaries where each dictionary contains a fitacf record (from pydarn.read_fitacf())
parameter: str
    The parameter to be detrended
    Default: 'both' (Velocity and SNR)
    Options: 'both', 'v', 'p_l'
window_length: int
    Length of the detrending low-pass filter in seconds
detrend_type: str
    Type of detrending to be used
    Default: 'mean' - subtract a running mean of window_length
    Option: 'savgol' - subtract a Savitsky-Golay filter instead
**kwargs:
    Optional inputs for detrending using scipy's `savgol_filter()`
    E.g., polyorder, mode
    Defaults are polyorder = 2, mode = 'interp' (edge truncation)

Returns
-------
fitacf_data_detrended: List[dict]
    Copy of input dmap data with detrended data substituted
"""

Approval

Number of approvals: 2

Test

Here is an example with SAS widebeam observations of Pc5 ULF waves on 2025-12-16:

import pydarn
import matplotlib.pyplot as plt
from pydarn.utils.detrend import Detrend

fitacf = '20251216.0000.27.sas.a.despecked.fitacf'
dmap_og, _ = pydarn.read_fitacf(fitacf)

dmap_detrended = Detrend.detrend_fitacf(dmap_og, parameter='both', window_length=600)

fig, ax = plt.subplots(4, 1)
pydarn.RTP.plot_range_time(dmap_og, beam_num=7, groundscatter=True, parameter='v',
                           zmax=500, zmin=-500, date_fmt='%H:%M', ax=ax[0],
                           range_estimation=pydarn.RangeEstimation.RANGE_GATE,
                           colorbar_label='LOS velocity\n[m/s]')

pydarn.RTP.plot_range_time(dmap_detrended, beam_num=7, groundscatter=True, parameter='v',
                           zmax=500, zmin=-500, date_fmt='%H:%M', ax=ax[1],
                           range_estimation=pydarn.RangeEstimation.RANGE_GATE,
                           colorbar_label='LOS velocity\n[m/s]')

pydarn.RTP.plot_range_time(dmap_og, beam_num=7, groundscatter=True, parameter='p_l',
                           zmax=0, zmin=50, date_fmt='%H:%M', ax=ax[2],
                           range_estimation=pydarn.RangeEstimation.RANGE_GATE,
                           colorbar_label='SNR [dB]')

pydarn.RTP.plot_range_time(dmap_detrended, beam_num=7, groundscatter=True, parameter='p_l',
                           zmax=0, zmin=50, date_fmt='%H:%M', ax=ax[3],
                           range_estimation=pydarn.RangeEstimation.RANGE_GATE,
                           colorbar_label='SNR [dB]')
[ax[i].set_xticklabels([]) for i in range(len(ax)-1)]
fig.autofmt_xdate()
ax[0].set_title('Original')
ax[1].set_title('Detrended')
ax[2].set_title('Original')
ax[3].set_title('Detrended')
fig.tight_layout()
plt.show()
image

Same code but with the Savitsky-Golay low-pass:

dmap_detrended = Detrend.detrend_fitacf(dmap_og, parameter='both', window_length=600, detrend_type='savgol')
image

@KhanalKrishna
Copy link
Contributor

Thanks, Daniel.
It works.
Test below:

Matplotlib version: 3.10.8

import pydarn
import matplotlib.pyplot as plt
from pydarn.utils.detrend import Detrend

fitacf = '20130202.2001.00.inv.fitacf'
dmap_og, _ = pydarn.read_fitacf(fitacf)

dmap_detrended = Detrend.detrend_fitacf(dmap_og, parameter='both', window_length=600)

fig, ax = plt.subplots(4, 1)
pydarn.RTP.plot_range_time(dmap_og, beam_num=7, groundscatter=True, parameter='v',ymin= 0, ymax= 50,
zmax=500, zmin=-500, date_fmt='%H:%M', ax=ax[0],
range_estimation=pydarn.RangeEstimation.RANGE_GATE,
colorbar_label='LOS velocity\n[m/s]')

pydarn.RTP.plot_range_time(dmap_detrended, beam_num=7, groundscatter=True, parameter='v',ymin= 0, ymax= 50,
zmax=500, zmin=-500, date_fmt='%H:%M', ax=ax[1],
range_estimation=pydarn.RangeEstimation.RANGE_GATE,
colorbar_label='LOS velocity\n[m/s]')

pydarn.RTP.plot_range_time(dmap_og, beam_num=7, groundscatter=True, parameter='p_l',ymin= 0, ymax= 50,
zmax=0, zmin=50, date_fmt='%H:%M', ax=ax[2],
range_estimation=pydarn.RangeEstimation.RANGE_GATE,
colorbar_label='SNR [dB]')

pydarn.RTP.plot_range_time(dmap_detrended, beam_num=7, groundscatter=True, parameter='p_l',ymin= 0, ymax= 50,
zmax=0, zmin=50, date_fmt='%H:%M', ax=ax[3],
range_estimation=pydarn.RangeEstimation.RANGE_GATE,
colorbar_label='SNR [dB]')
[ax[i].set_xticklabels([]) for i in range(len(ax)-1)]
fig.autofmt_xdate()
ax[0].set_title('Original')
ax[1].set_title('Detrended')
ax[2].set_title('Original')
ax[3].set_title('Detrended')
fig.tight_layout()
plt.show()

image

With dtrend_type= 'savgol';
dmap_detrended = Detrend.detrend_fitacf(dmap_og, parameter='both', window_length=600, dtrend_type= 'savgol')
image

@carleyjmartin
Copy link
Collaborator

Why is the detrended SNR so much lower in magnitude?

@billetd
Copy link
Contributor Author

billetd commented Feb 12, 2026

Why is the detrended SNR so much lower in magnitude?

Detrending removes background fluctuations of the timescale you choose for the filter, highlighting only fluctuations shorter. The SNR going to zero just means that it's varying slowly enough to get filtered out. It's meant to highlight periodicity, but could also be used to identify transient spikes.

Copy link
Collaborator

@carleyjmartin carleyjmartin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Testing was fine, no issues.

You will want to import this feature into the init.py file.
from .utils.detrend import Detrend
Then you don't have to type pydarn.utils.detrend.Detrend.detrend_fitacf() when using this.

I'm still not sure why this can't just go in the filters module still anyway. You can have another class Detrend in there if you want, I just don't want us to blow up with modules with one thing in each of them.

Also, this needs adding with examples to documentation pls.

@@ -0,0 +1,268 @@
import copy
import pydarn
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're gonna want to do this as
from pydarn import SuperDARNRadars

I think that's all you want, so instead of importing the entire package again, just add this below.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants