liebharc / basic_dsp

Basic DSP vector operations for Rust.
Apache License 2.0
43 stars 5 forks source link

Convolution #50

Open HerrMuellerluedenscheid opened 1 year ago

HerrMuellerluedenscheid commented 1 year ago

Hi,

I have the following test main.rs which which I tried to convole the raised sin with some signal. When this is carried out, I see that apparently only the first three floats slightly changed somewhere in the lower digits. I was assuming that convolve behaves somewhat like numpys convolve but I guess I am mistaken here. How can I convolve the entire dsp_vec with the function in the example below?

Looking forward hearing from you

use basic_dsp::{SingleBuffer, ToRealVector};
use basic_dsp::conv_types::*;
use basic_dsp::*;
use basic_dsp::numbers::{Float, RealNumber};
use rand::prelude::*;

const INIT_VAL_RANGE: (f64, f64) = (-100.0, 100.0);

fn create_pseudo_random_data(data_set_size: usize, seed: usize) -> Vec<f64> {
    let mut vec: Vec<f64> = vec![0.0; data_set_size];
    let mut rng: StdRng = SeedableRng::seed_from_u64(seed as u64);
    for n in &mut vec {
        *n = rng.gen_range(std::ops::Range { start: INIT_VAL_RANGE.0, end: INIT_VAL_RANGE.1 });
    }
    vec
}

fn main() {
    let data_set_size = 200;
    let mut dsp_vec = create_pseudo_random_data(data_set_size, 5798734).to_real_time_vec();
    for i in 0..data_set_size / 2 {
        dsp_vec.data[i] = (i as f64 / 10.).sin() * 5.0;
    }
    println!("{:?}", &dsp_vec.data[1..10]);

    let mut buffer = SingleBuffer::new();
    let function = RaisedCosineFunction::new(0.25);
    dsp_vec.convolve(&mut buffer, &function, 2.0, 12);
    let vec: Vec<f64> = dsp_vec.into();
    println!("{:?}", &vec[1..10]);
}
liebharc commented 1 year ago

Possibly the ratio factor isn't giving you the impulse response you need. If you print the values which you get in this configuration then you see that only the center value is 1 and all the other values are close to 0. That's a dirac delta, which makes the convolution a noop (that is: the output equals the input). I added a little bit of code to print a table to illustrate the function values:

use basic_dsp::{SingleBuffer, ToRealVector};
use basic_dsp::conv_types::*;
use basic_dsp::*;
use basic_dsp::numbers::{Float, RealNumber};
use rand::prelude::*;

const INIT_VAL_RANGE: (f64, f64) = (-100.0, 100.0);

fn create_pseudo_random_data(data_set_size: usize, seed: usize) -> Vec<f64> {
    let mut vec: Vec<f64> = vec![0.0; data_set_size];
    let mut rng: StdRng = SeedableRng::seed_from_u64(seed as u64);
    for n in &mut vec {
        *n = rng.gen_range(std::ops::Range { start: INIT_VAL_RANGE.0, end: INIT_VAL_RANGE.1 });
    }
    vec
}

fn main() {

    let mut buffer = SingleBuffer::new();
    let function = RaisedCosineFunction::new(0.25);

    let len = 12;
    let ratio = 2.0;
    let mut j: f64 = -(len as f64);
    while j <= len as f64 {
        let value = basic_dsp::conv_types::RealImpulseResponse::calc(&function, j * ratio);
        println!("{:?}, {:?}", j, value);
        j = j + 1.0;
    }

    let data_set_size = 200;
    let mut dsp_vec = create_pseudo_random_data(data_set_size, 5798734).to_real_time_vec();
    for i in 0..data_set_size / 2 {
        dsp_vec.data[i] = (i as f64 / 10.).sin() * 5.0;
    }
    let org = dsp_vec.clone();
    let range = 0..10;
    println!("{:?}", &dsp_vec.data[range.clone()]);

    dsp_vec.convolve(&mut buffer, &function as &dyn RealImpulseResponse<f64>, ratio, len);
    println!("{:?}", &dsp_vec[range.clone()]);
    dsp_vec.sub(&org);
    println!("{:?}", &dsp_vec[range.clone()]);
}

which then gives this:

-12.0, 2.7259942884750876e-19  
-11.0, -2.894725366663278e-33  
-10.0, -3.9375473055751276e-19 
-9.0, 2.6853095571869394e-34   
-8.0, 6.187574337332342e-19    
-7.0, -3.480956833390477e-34   
-6.0, -1.1137633807198216e-18  
-5.0, 4.97279547627211e-34     
-4.0, 2.5987812216795837e-18   
-3.0, -8.951031857289796e-34   
-2.0, -1.2993906108397918e-17  
-1.0, -3.061616997868383e-17   
0.0, 1.0                       
1.0, -3.061616997868383e-17    
2.0, -1.2993906108397918e-17   
3.0, -8.951031857289796e-34    
4.0, 2.5987812216795837e-18    
5.0, 4.97279547627211e-34      
6.0, -1.1137633807198216e-18   
7.0, -3.480956833390477e-34    
8.0, 6.187574337332342e-19     
9.0, 2.6853095571869394e-34    
10.0, -3.9375473055751276e-19  
11.0, -2.894725366663278e-33   
12.0, 2.7259942884750876e-19   

I think the documentation is lacking in this regard and could use some improvement. Also, the usability of the API and the test coverage for vectors in real number space could be better. Unfortunately, I don't have the time and energy to continue to maintain this repo on my own. I would mark this issue as "help wanted", but don't expect that anyone will pick it up - at least not anytime soon. So, apologies if this response isn't as helpful as you might have hoped it to be.

Viele Grüße nach Berlin

liebharc commented 1 year ago

If someone wants to work on this issue:

HerrMuellerluedenscheid commented 1 year ago

❤️ Besten Dank! That was very exhaustive. I'll try to dedicate some time in expanding on the documentation and create a little example once I have achieved what I wanted.