6.10.94-stable
In this page
  1. Generate the time series
  2. Quick access to FFT model
  3. Quick Extrapolate method

FFT (Fast Fourier Transform)

FFT is a mathematical tool used to decompose a signal to its underlying base frequencies (known as harmonics). This is very helpful to detect the periodicity aspect in a signal.

To use FFT in greycat, you need to add the algebra library and use transforms module as following:

@library("algebra");
use transforms;

Generate the time series

In order to show a usage example of FFT, let’s first generate a time series at random and store it in a nodeTime as following:

use util;

fn main() {
    var ts = nodeTime::new();
    //Random number generator with a fixed seed
    var rand = Random::new();
    rand.setSeed(1234);

    var starting = time::now().to(DurationUnit::microseconds);
    var avgT = 10000; //generate a time point every 10 ms
    var N = 17;

    var offset = 0;
    println("Displaying original time series");
    for (var i = 0; i < N; i++) {
        offset = offset + avgT;
        var v = rand.uniformf(0.0,100.0);
        ts.setAt(time::new(starting + offset, DurationUnit::microseconds), v);
        println("${i} - ${offset} => ${v}");
    }
    println("");
}

Quick access to FFT model

If you don’t want to pilot the FFT by yourself, and you care about getting a quick FFT model of a time series, we have the easy to use following methods:

var fft_model = FFTModel::train(ts, ts.firstTime()!!, ts.lastTime()!!);

This creates automatically an FFT decomposition of a time series, you just need to provide the time series, first time-point where to start the analysis and last time-point. Inside the model you will find the previously described best_size, frequency_complex, frequency_table etc

Quick Extrapolate method

The previously generated fft_model can be used to generate a table of time series as following:

var result = fft_model.extrapolate(start_time, end_time, low_pass_filter_ratio, nb_harmonics, skip_elements);
  • low_pass_filter_ratio: non mandatory attribute, a ratio from 0 to 1.0 (keep all frequencies), low_pass_filter_ratio=0.95 keeps 95% of the lowest frequencies
  • nb_harmonics: non mandatory attribute, select the max number of frequencies to take into account
  • skip_elements non mandatory attributes, if skip_elements=0 no elements in the generation are skipped, if you specify the skip elements number the generated table will be smaller in size (not all elements are generated).

Manual FFT usage

Preparing the attributes

If you want more advanced usage of FFT, you need to do the following steps:

    var ts_start = ts.firstTime()!!; //Get the first time point of the time series
    var ts_end = ts.lastTime()!!;   //Get the last time point of the time series
    var size = ts.rangeSize(ts_start, ts_end); //get the total number of elements in the time series

Now we have the start time, end time and the size of the time series. FFT algorithm is optimized to run better on specific sizes, so to get the optimize size run the following:

    var best_size = FFT::get_next_fast_size(size); //gets an optimized size for FFT
    var sampling_step = (ts_end - ts_start).to(DurationUnit::microseconds) / (best_size - 1);
    var sampling_duration = duration::new(sampling_step, DurationUnit::microseconds);
    info("time series has ${size} time points => next best size is: ${best_size}, sampling duration is: ${sampling_step} us");

To create a forward transform FFT, we need to create the following:

    var fft = FFT::new(best_size,false); //Creates a forward transform fourier
    var time_complex = Tensor {};   //Create a time complex tensor to store the time series

Note that the FFT object is un-serializable, so it can’t be stored to disk. It contains however an array of size N of precalculated cosines (to accelerate the computation). So if you need to run fourier transform over windows of same length, think to re-use this object for performance.

Running forward transform

To fill the time_complex tensor, we need to iterate on the nodeTime at a fixed sampling rate, we do so using the following:

    //Initialize the time_complex with c128 tensor type, and dimension best_size
    time_complex.init(TensorType::c128, [best_size]); 

    //Initialize a nodeTimeCursor to iterate on the nodeTime
    var nc: nodeTimeCursor = nodeTimeCursor::new();
    nc.init(ts);

    println("Displaying sampled time series");
    for (var i = 0; i < best_size; i++) {
        //set the time_complex tensor 
        time_complex.set([i], nc.current() as float);

        println("${i} - ${nc.currentTime()!!.to(DurationUnit::microseconds) - starting} => ${nc.current()}");
        nc.skip_duration(sampling_duration);
    }
    println("");

Last step, in order to convert time_complex tensor to frequency_complex we do the following:

    //Convert time complex to frequency complex
    var frequency_complex = Tensor {};
    fft.transform(time_complex, frequency_complex);

After the fft.transform call you get the frequency tensor as a complex object containing the FFT decomposition of the time_complex tensor.

Frequency table

FFT class contains several helpers to improve the usage of the fft_complex tensor in gcl. One such helper is the get_frequency_table:

    var frequency_table = FFT::get_frequency_table(frequency_complex,sampling_duration);

This method converts the frequency_complex tensor into a table useful for display in the frontend. The table contains the following columns:

  • Column 0: Frequency in hz: the frequency of the harmonic in Hz, useful to be displayed in X axis
  • Column 1: Frequency abs power: the absolute power of the harmonic, useful to be displayed in Y axis
  • Column 2: phase angle (in radians): the phase in radians of the harmonic
  • Column 3: position: an indexed position of the frequency in the spectrum, useful to keep track of the position if the table gets sorted
  • Column 4: period (duration): The inverse of the frequency in Hz (column 0), it’s returned as a greycat duration type
  • Column 5: power in db: the absolute power of the harmonic calculated in decibel, 20log10(column1)
  • Column 6: cumulated ratio: the cumulated ratio of the harmonics so far, at the last frequency this should be 1.0. The cumulated ratio helps to decide where to cut the frequencies in a low_pass_filter to keep let’s say 95% of the signal.

Low pass filter

    var lowpass = FFT::get_low_pass_filter_size(frequency_complex, 0.95);
    println(lowpass);

This method gives the index of where to cut the frequency spectrum to keep 95% of the signal. In signal processing it’s called low pass filtering.