I have a continuous value for which I&#39;d like to calculate an [exponential moving average](http://en.wikipedia.org/wiki/Exponential_moving_average#Exponential_moving_average). Normally I&#39;d just use the standard formula for this:

* S&lt;sub&gt;n&lt;/sub&gt; = &amp;alpha;Y + (1-&amp;alpha;)S&lt;sub&gt;n-1&lt;/sub&gt;

where S&lt;sub&gt;n&lt;/sub&gt; is the new average, &amp;alpha; is the alpha, Y is the sample, and S&lt;sub&gt;n-1&lt;/sub&gt; is the previous average.

Unfortunately, due to various issues I don&#39;t have a consistent sample time. I may know I can sample at the most, say, once per millisecond, but due to factors out of my control, I may not be able to take a sample for several milliseconds at a time. A likely more common case, however, is that I simple sample a bit early or late: instead of sampling at 0, 1 and 2 ms. I sample at 0, 0.9 and 2.1 ms. I do anticipate that, regardless of delays, my sampling frequency will be far, far above the Nyquist limit, and thus I need not worry about aliasing.

I reckon that I can deal with this in a more-or-less reasonable way by varying the alpha appropriately, based on the length of time since the last sample.

Part of my reasoning that this will work is that the EMA &quot;interpolates linearly&quot; between the previous data point and the current one. If we consider calculating an EMA of the following list of samples at intervals t: [0,1,2,3,4]. We should get the same result if we use interval 2t, where the inputs become [0,2,4], right? If the EMA had assumed that, at t&lt;sub&gt;2&lt;/sub&gt; the value had been 2 since t&lt;sub&gt;0&lt;/sub&gt;, that would be the same as the interval t calculation calculating on [0,2,2,4,4], which it&#39;s not doing. Or does that make sense at all?

Can someone tell me how to vary the alpha appropriately? &quot;Please show your work.&quot; I.e., show me the math that proves that your method really is doing the right thing.

Exponential Moving Average Sampled at Varying Times

Lets assume we have a dataset which might be given approximately by

    import numpy as np
    x = np.linspace(0,2*np.pi,100)
    y = np.sin(x) + np.random.random(100) * 0.2

Therefore we have a variation of 20% of the dataset. My first idea was to use the UnivariateSpline function of scipy, but the problem is that this does not consider the small noise in a good way. If you consider the frequencies, the background is much smaller than the signal, so a spline only of the cutoff might be an idea, but that would involve a back and forth fourier transformation, which might result in bad behaviour.
Another way would be a moving average, but this would also need the right choice of the delay.

Any hints/ books or links how to tackle this problem?

![example][1]
    


  [1]: http://i.stack.imgur.com/bcaCP.png

How to smooth a curve in the right way?

Can anyone please clearly explain the difference between 1D, 2D, and 3D convolutions in convolutional neural networks (in deep learning) with the use of examples?

Intuitive understanding of 1D, 2D, and 3D convolutions in convolutional neural networks

I have an FFT result. These are stored in two `double` arrays: a real part array and an imaginary part array. How do I determine the frequencies that correspond to each element in these arrays?

In other words, I would like have create an array that stores the frequencies for each real and imaginary component of my FFT.


How do I obtain the frequencies of each value in an FFT?

I am trying to create a graphical spectrum analyzer in python.

I am currently reading 1024 bytes of a 16 bit dual channel 44,100 Hz sample rate audio stream and averaging the amplitude of the 2 channels together.  So now I have an array of 256 signed shorts.  I now want to preform a fft on that array, using a module like numpy, and use the result to create the graphical spectrum analyzer, which, to start will just be 32 bars.

I have read the wikipedia articles on Fast Fourier Transform and Discrete Fourier Transform but I am still unclear of what the resulting array represents.  This is what the array looks like after I preform an fft on my array using numpy:

 

       [ -3.37260500e+05 +0.00000000e+00j   7.11787022e+05 +1.70667403e+04j
       4.10040193e+05 +3.28653370e+05j   9.90933073e+04 +1.60555003e+05j
       2.28787050e+05 +3.24141951e+05j   2.09781047e+04 +2.31063376e+05j
      -2.15941453e+05 +1.63773851e+05j  -7.07833051e+04 +1.52467334e+05j
      -1.37440802e+05 +6.28107674e+04j  -7.07536614e+03 +5.55634993e+03j
      -4.31009964e+04 -1.74891657e+05j   1.39384348e+05 +1.95956947e+04j
       1.73613033e+05 +1.16883207e+05j   1.15610357e+05 -2.62619884e+04j
      -2.05469722e+05 +1.71343186e+05j  -1.56779748e+04 +1.51258101e+05j
      -2.08639913e+05 +6.07372799e+04j  -2.90623668e+05 -2.79550838e+05j
      -1.68112214e+05 +4.47877871e+04j  -1.21289916e+03 +1.18397979e+05j
      -1.55779104e+05 +5.06852464e+04j   1.95309737e+05 +1.93876325e+04j
      -2.80400414e+05 +6.90079265e+04j   1.25892113e+04 -1.39293422e+05j
       3.10709174e+04 -1.35248953e+05j   1.31003438e+05 +1.90799303e+05j...

I am wondering what exactly these numbers represent and how I would convert these numbers into a percentage of a height for each of the 32 bars.  Also, should I be averaging the 2 channels together?

Analyze audio using Fast Fourier Transform

I need some help understanding the output of the DFT/FFT computation.

I&#39;m an experienced software engineer and need to interpret some smartphone accelerometer readings, such as finding the principal frequencies. Unfortunately, I slept through most of my college EE classes fifteen years ago, but I&#39;ve been reading up on DFT and FFT for the last several days (to little avail, apparently). 

Please, no responses of &quot;go take an EE class&quot;. I&#39;m actually planning to do that if my employer will pay me. :)

So here is my problem:


I&#39;ve captured a signal at 32 Hz. Here is a 1 second sample of 32 points, which I&#39;ve charted in Excel.

![enter image description here][1]


I then got some [FFT code][2] written in Java from Columbia University (after following the suggestions in a post on &quot;[Reliable and fast FFT in Java][3]&quot;).

The output of this program is as follows. I believe it is running an in-place FFT, so it re-uses the same buffer for both input and output.

    Before: 
    
    Re: [0.887  1.645  2.005  1.069  1.069  0.69  1.046  1.847  0.808  0.617  0.792  1.384  1.782  0.925  0.751  0.858  0.915  1.006  0.985  0.97  1.075  1.183  1.408  1.575  1.556  1.282  1.06  1.061  1.283  1.701  1.101  0.702  ]
    
    Im: [0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  ]
    
    After: 
    
    Re: [37.054  1.774  -1.075  1.451  -0.653  -0.253  -1.686  -3.602  0.226  0.374  -0.194  -0.312  -1.432  0.429  0.709  -0.085  0.0090  -0.085  0.709  0.429  -1.432  -0.312  -0.194  0.374  0.226  -3.602  -1.686  -0.253  -0.653  1.451  -1.075  1.774  ]
    
    Im: [0.0  1.474  -0.238  -2.026  -0.22  -0.24  -5.009  -1.398  0.416  -1.251  -0.708  -0.713  0.851  1.882  0.379  0.021  0.0  -0.021  -0.379  -1.882  -0.851  0.713  0.708  1.251  -0.416  1.398  5.009  0.24  0.22  2.026  0.238  -1.474  ]


So, at this point, I can&#39;t make heads or tails of the output. I understand the DFT concepts, such as the real portion being the amplitudes of the component cosine waves and the imaginary portion being the amplitudes of the component sine waves. I can also follow this diagram from the great book &quot;[The Scientist and Engineer&#39;s Guide to Digital Signal Processing][4]&quot;:
![enter image description here][5]


So my specific questions are:

1. From the output of the FFT, how do I find the &quot;most occurring frequencies&quot;? This is part of my analysis of my accelerometer data. Should I read the real (cosine) or imaginary (sine) arrays?

2. I have a 32-point input in the time domain. Shouldn&#39;t the output of the FFT be a 16-element array for reals and a 16-element array for imaginary? Why does the program give me real and imaginary array outputs both of size 32?

3. Related to the previous question, how do I parse the indexes in the output arrays? Given my input of 32 samples sampled at 32 Hz, my understanding is that a 16-element array output should have its index uniformly spread up to 1/2 the sampling rate (of 32 Hz), so am I correct in understanding that each element of the array represents (32 Hz * 1/2) / 16 = 1 Hz?

4. Why does the FFT output have negative values? I thought the values represent amplitudes of a sinusoid. For example, the output of Real[ 3 ] = -1.075 should mean an amplitude of -1.075 for a cosine wave of frequency 3. Is that right? How can an amplitude be negative?


  [1]: http://i.stack.imgur.com/IHj1I.jpg
  [2]: http://www.ee.columbia.edu/~ronw/code/MEAPsoft/doc/html/FFT_8java-source.html
  [3]: https://stackoverflow.com/questions/3287518/reliable-and-fast-fft-in-java
  [4]: http://www.dspguide.com/pdfbook.htm
  [5]: http://i.stack.imgur.com/MDI8J.jpg

Understanding FFT output

UPDATE:

I found a Scipy Recipe based in this question! So, for anyone interested, go straight to: [Contents &#187; Signal processing &#187; Butterworth Bandpass][1]

----------

I&#39;m having a hard time to achieve what seemed initially a simple task of implementing a Butterworth band-pass filter for 1-D numpy array (time-series).

The parameters I have to include are the sample_rate, cutoff frequencies IN HERTZ and possibly order (other parameters, like attenuation, natural frequency, etc. are more obscure to me, so any &quot;default&quot; value would do).

What I have now is this, which seems to work as a high-pass filter but I&#39;m no way sure if I&#39;m doing it right:

    def butter_highpass(interval, sampling_rate, cutoff, order=5):
        nyq = sampling_rate * 0.5
        
        stopfreq = float(cutoff)
        cornerfreq = 0.4 * stopfreq  # (?)
        
        ws = cornerfreq/nyq
        wp = stopfreq/nyq
        
        # for bandpass:
        # wp = [0.2, 0.5], ws = [0.1, 0.6]
        
        N, wn = scipy.signal.buttord(wp, ws, 3, 16)   # (?)
        
        # for hardcoded order:
        # N = order
        
        b, a = scipy.signal.butter(N, wn, btype=&#39;high&#39;)   # should &#39;high&#39; be here for bandpass?
        sf = scipy.signal.lfilter(b, a, interval)
        return sf

![enter image description here][2]

The docs and examples are confusing and obscure, but I&#39;d like to implement the form presented in the commend marked as &quot;for bandpass&quot;. The question marks in the comments show where I just copy-pasted some example without understanding what is happening.

I am no electrical engineering or scientist, just a medical equipment designer needing to perform some rather straightforward bandpass filtering on EMG signals.

  [1]: http://scipy-cookbook.readthedocs.io/items/ButterworthBandpass.html
  [2]: http://i.stack.imgur.com/CtpRF.png

How to implement band-pass Butterworth filter with Scipy.signal.butter

I am trying to filter a noisy heart rate signal with python. Because heart rates should never be above about 220 beats per minute, I want to filter out all noise above 220 bpm. I converted 220/minute into 3.66666666 Hertz and then converted that Hertz to rad/s to get 23.0383461 rad/sec.

The sampling frequency of the chip that takes data is 30Hz so I converted that to rad/s to get 188.495559 rad/s.

After looking up some stuff online I found some functions for a bandpass filter that I wanted to make into a lowpass. [Here is the link the bandpass code][1], so I converted it to be this:

    from scipy.signal import butter, lfilter
    from scipy.signal import freqs
    
    def butter_lowpass(cutOff, fs, order=5):
        nyq = 0.5 * fs
        normalCutoff = cutOff / nyq
        b, a = butter(order, normalCutoff, btype=&#39;low&#39;, analog = True)
        return b, a
    
    def butter_lowpass_filter(data, cutOff, fs, order=4):
        b, a = butter_lowpass(cutOff, fs, order=order)
        y = lfilter(b, a, data)
        return y
    
    cutOff = 23.1 #cutoff frequency in rad/s
    fs = 188.495559 #sampling frequency in rad/s
    order = 20 #order of filter
    
    #print sticker_data.ps1_dxdt2
    
    y = butter_lowpass_filter(data, cutOff, fs, order)
    plt.plot(y)


I am very confused by this though because I am pretty sure the butter function takes in the cutoff and sampling frequency in rad/s but I seem to be getting a weird output. Is it actually in Hz?

Secondly, what is the purpose of these two lines:

        nyq = 0.5 * fs
        normalCutoff = cutOff / nyq

I know it&#39;s something about normalization but I thought the nyquist was 2 times the sampling requency, not one half. And why are you using the nyquist as a normalizer?

Can someone explain more about how to create filters with these functions?

I plotted the filter using:

    w, h = signal.freqs(b, a)
    plt.plot(w, 20 * np.log10(abs(h)))
    plt.xscale(&#39;log&#39;)
    plt.title(&#39;Butterworth filter frequency response&#39;)
    plt.xlabel(&#39;Frequency [radians / second]&#39;)
    plt.ylabel(&#39;Amplitude [dB]&#39;)
    plt.margins(0, 0.1)
    plt.grid(which=&#39;both&#39;, axis=&#39;both&#39;)
    plt.axvline(100, color=&#39;green&#39;) # cutoff frequency
    plt.show()

and got this which clearly does not cut-off at 23 rad/s:

![result][2]

  [1]: http://wiki.scipy.org/Cookbook/ButterworthBandpass
  [2]: http://i.imgur.com/weqB1sw.png

Creating lowpass filter in SciPy - understanding methods and units

All the FFT implementations we have come across result in complex values (with real and imaginary parts), even if the input to the algorithm was a discrete set of real numbers (integers).

Is it not possible to represent frequency domain in terms of real numbers only?

Why does FFT produce complex numbers instead of real numbers?

I&#39;m trying to read an image from an electrocardiography and detect each one of the main waves in it (P wave, QRS complex and T wave). I can read the image and get a vector (like `(4.2; 4.4; 4.9; 4.7; ...)`). I need an algorithm that can walk through this vector and detect when each of these waves start and end. An example:

![alt text][1]

Would be easy if they always had the same size, or if I knew how many waves the ECG has in advance. Given the wave:

![alt text][2]

I extract the vector:

    [0; 0; 20; 20; 20; 19; 18; 17; 17; 17; 17; 17; 16; 16; 16; 16; 16; 16; 16; 17; 17; 18; 19; 20; 21; 22; 23; 23; 23; 25; 25; 23; 22; 20; 19; 17; 16; 16; 14; 13; 14; 13; 13; 12; 12; 12; 12; 12; 11; 11; 10; 12; 16; 22; 31; 38; 45; 51; 47; 41; 33; 26; 21; 17; 17; 16; 16; 15; 16; 17; 17; 18; 18; 17; 18; 18; 18; 18; 18; 18; 18; 17; 17; 18; 19; 18; 18; 19; 19; 19; 19; 20; 20; 19; 20; 22; 24; 24; 25; 26; 27; 28; 29; 30; 31; 31; 31; 32; 32; 32; 31; 29; 28; 26; 24; 22; 20; 20; 19; 18; 18; 17; 17; 16; 16; 15; 15; 16; 15; 15; 15; 15; 15; 15; 15; 15; 15; 14; 15; 16; 16; 16; 16; 16; 16; 16; 16; 16; 15; 16; 15; 15; 15; 16; 16; 16; 16; 16; 16; 16; 16; 15; 16; 16; 16; 16; 16; 15; 15; 15; 15; 15; 16; 16; 17; 18; 18; 19; 19; 19; 20; 21; 22; 22; 22; 22; 21; 20; 18; 17; 17; 15; 15; 14; 14; 13; 13; 14; 13; 13; 13; 12; 12; 12; 12; 13; 18; 23; 30; 38; 47; 51; 44; 39; 31; 24; 18; 16; 15; 15; 15; 15; 15; 15; 16; 16; 16; 17; 16; 16; 17; 17; 16; 17; 17; 17; 17; 18; 18; 18; 18; 19; 19; 20; 20; 20; 20; 21; 22; 22; 24; 25; 26; 27; 28; 29; 30; 31; 32; 33; 32; 33; 33; 33; 32; 30; 28; 26; 24; 23; 23; 22; 20; 19; 19; 18; 17; 17; 18; 17; 18; 18; 17; 18; 17; 18; 18; 17; 17; 17; 17; 16; 17; 17; 17; 18; 18; 17; 17; 18; 18; 18; 19; 18; 18; 17; 18; 18; 17; 17; 17; 17; 17; 18; 17; 17; 18; 17; 17; 17; 17; 17; 17; 17; 18; 17; 17; 18; 18; 18; 20; 20; 21; 21; 22; 23; 24; 23; 23; 21; 21; 20; 18; 18; 17; 16; 14; 13; 13; 13; 13; 13; 13; 13; 13; 13; 12; 12; 12; 16; 19; 28; 36; 47; 51; 46; 40; 32; 24; 20; 18; 16; 16; 16; 16; 15; 16; 16; 16; 17; 17; 17; 18; 17; 17; 18; 18; 18; 18; 19; 18; 18; 19; 20; 20; 20; 20; 20; 21; 21; 22; 22; 23; 25; 26; 27; 29; 29; 30; 31; 32; 33; 33; 33; 34; 35; 35; 35; 0; 0; 0; 0;]

I would like to detect, for example:

 - P wave in `[19 - 37]`.
 - QRS complex in `[51 - 64]`.
 - etc.

  [1]: http://i.stack.imgur.com/U4t3t.jpg
  [2]: http://i.stack.imgur.com/BHOaG.jpg

How to detect patterns in (electrocardiography) waves?

How can the tempo/BPM of a song be determined programmatically?  What algorithms are commonly used, and what considerations must be made?

How to detect the BPM of a song in php

Where can I find a free, very quick, and reliable implementation of FFT in C#?

That can be used in a product? Or are there any restrictions?


An implementation of the fast Fourier transform (FFT) in C#

Has anybody used the `Apple FFT` for an iPhone app yet or know where I might find a sample application as to how to use it?  I know that Apple has some sample code posted, but I&#39;m not really sure how to implement it into an actual project.

Using the Apple FFT and Accelerate Framework

---

**Update:** The best performing algorithm *so far* [**is this one**](https://stackoverflow.com/questions/22583391/peak-recognition-in-realtime-timeseries-data/22640362#22640362).

---

*This question explores robust algorithms for detecting sudden peaks in real-time timeseries data.*


Consider the following example data:

[![Plot of data][1]][1]

&lt;sub&gt;Example of this data is in Matlab format (but this question is not about the language but about the algorithm):&lt;/sub&gt;

    p = [1 1 1.1 1 0.9 1 1 1.1 1 0.9 1 1.1 1 1 0.9 1 1 1.1 1 1 1 1 1.1 0.9 1 1.1 1 1 0.9, ...
         1 1.1 1 1 1.1 1 0.8 0.9 1 1.2 0.9 1 1 1.1 1.2 1 1.5 1 3 2 5 3 2 1 1 1 0.9 1 1, ... 
         3 2.6 4 3 3.2 2 1 1 0.8 4 4 2 2.5 1 1 1];

You can clearly see that there are three large peaks and some small peaks. This dataset is a specific example of the class of timeseries datasets that the question is about. This class of datasets has two general features:

 1. There is basic noise with a general mean
 2. There are large &#39;*peaks*&#39; or &#39;*higher data points*&#39; that significantly deviate from the noise. 

Let&#39;s also assume the following:

 - The width of the peaks cannot be determined beforehand
 - The height of the peaks significantly deviates from the other values
 - The algorithm updates in realtime (so updates with each new datapoint)

For such a situation, a boundary value needs to be constructed which triggers signals. However, the boundary value cannot be static and must be determined realtime based on an algorithm.

---

**My Question: what is a good algorithm to calculate such thresholds in realtime?** Are there specific algorithms for such situations? What are the most well-known algorithms?

---

&lt;sub&gt;Robust algorithms or useful insights are all highly appreciated. (can answer in any language: it&#39;s about the algorithm)&lt;/sub&gt;


  [1]: https://i.stack.imgur.com/yUeKr.jpg

Peak signal detection in realtime timeseries data

I am testing the WebRTC AGC but I must be doing something wrong because the signal just passes through unmodified.

Here&#39;s how I create and initialize the AGC:

    agcConfig.compressionGaindB = 9;
    agcConfig.limiterEnable = 1;
    agcConfig.targetLevelDbfs = 9;   /* 9dB below full scale */

    WebRtcAgc_Create(&amp;agc);
    WebRtcAgc_Init(agc, minLevel, maxLevel, kAgcModeFixedDigital, 8000);
    WebRtcAgc_set_config(agc, agcConfig);

And then for each 10ms sample block I do the following:

    WebRtcAgc_Process(agc, micData, NULL, 80, micData, NULL, micLevelIn, &amp;micLevelOut, 0, &amp;saturationWarning);

Where `micLevelIn` is set to 0. 

Can somebody tell me what I&#39;m doing wrong? 

I expected that a full scale sine tone would be attenuated to the target DBFS level; and a low level sine tone (i.e. -30dBFS) would be amplified to match the target DBFS level. But that&#39;s not what I&#39;m seeing.

WebRTC AGC (Automatic Gain Control)

I want to calculate room noise level with the computer&#39;s microphone. I record noise as an audio file, but how can I calculate the noise dB level?

I don&#39;t know how to start!

All Signal Processing Solutions on Gang of Coders

Total of 16 Signal Processing Solutions