Noise estimation by polynomial approximation¶
The estimateSNR() function provides an algorithm
to estimate the SNR in a data set.
- PyAstronomy.pyasl.estimateSNR(x, y, xlen, deg=1, controlPlot=False, xlenMode='dataPoints')¶
Estimate Signal to Noise Ratio (SNR) in a data set.
This function provides a simple algorithm to estimate the SNR in a data set. The algorithm subdivides the data set into sections with a length determined by the xlen parameter. Each individual subsection is fitted using a polynomial of degree deg. The SNR is then computed by assuming that the resulting chi square value is properly distributed according to the chi-square distribution.
- Parameters
- xarray
The abscissa values.
- yarray
The ordinate values.
- xlenint or float
Length of the data subsection considered in the computation of the SNR. Whether xlen refers to the number of data points or a fixed subsection of the abscissa is determined by the xlenMode flag.
- degint, optional
Degree of the polynomial used to fit the subsection of the data set. The default is one.
- controlPlotboolean, optional
If True, a control plot will be shown to verify the validity of the estimate.
- xlenModestring, {“dataPoints”, “excerpt”, “all”}, optional
Determines whether xlen refers to data points or a fixed length scale on the abscissa. If ‘all’ is specified, all available data will be used in the estimation.
- Returns
- SNRdictionary
Contains the SNR estimate for the individual subsections (key ‘snrs’) and the final SNR estimate (mean of all ‘sub-SNRs’, key ‘SNR-Estimate’).
Example—Estimating SNR ratio¶
from __future__ import print_function, division
from PyAstronomy import pyasl
import numpy as np
# Number of data points
N = 10000
# Signal to noise ratio
SNR = 50.0
# Create some data with noise and a sinusoidal
# variation.
x = np.arange(N)
y = np.random.normal(0.0, 1.0/SNR, N) + 1.0
y += np.sin(x/500.0*2.*np.pi)*0.1
# Estimate the signal to noise ratio. Check whether the
# estimate fits the input...
# Use a chunk length of 20 data points, a polynomial of degree
# one, and produce a "control plot".
snrEsti = pyasl.estimateSNR(x, y, 20, deg=1, controlPlot=True)
print("Estimate of the SNR: ", snrEsti["SNR-Estimate"])
# Use a chunks with a length of 27, a polynomial of degree
# two, and produce a "control plot".
snrEsti = pyasl.estimateSNR(x, y, 27, deg=2, controlPlot=False, xlenMode="excerpt")
print("Estimate of the SNR: ", snrEsti["SNR-Estimate"])
