spectral_analysis/
spectral_analysis.rs

1//! Exemplo de análise espectral usando avila-fft
2//!
3//! Demonstra:
4//! - Geração de sinal composto (múltiplas frequências)
5//! - Aplicação de janelamento
6//! - FFT e análise do espectro
7//! - Detecção de picos de frequência
8
9use avila_fft::{Complex, FftPlanner, window};
10use std::f64::consts::PI;
11
12fn main() {
13    println!("=== Análise Espectral com Avila FFT ===\n");
14
15    // Parâmetros do sinal
16    let sample_rate = 1000.0; // Hz
17    let duration = 1.0; // segundos
18    let n = (sample_rate * duration) as usize; // 1000 amostras
19
20    // Garante potência de 2
21    let n = n.next_power_of_two();
22    println!("Amostras: {}", n);
23    println!("Taxa de amostragem: {} Hz", sample_rate);
24    println!("Resolução: {:.2} Hz\n", sample_rate / n as f64);
25
26    // Gera sinal composto: 50Hz + 120Hz + 300Hz + ruído
27    println!("Gerando sinal com 3 componentes:");
28    println!("  - 50 Hz (amplitude 1.0)");
29    println!("  - 120 Hz (amplitude 0.7)");
30    println!("  - 300 Hz (amplitude 0.3)");
31    println!("  + ruído branco\n");
32
33    let signal: Vec<f64> = (0..n)
34        .map(|i| {
35            let t = i as f64 / sample_rate;
36
37            // Componentes senoidais
38            let sig1 = 1.0 * (2.0 * PI * 50.0 * t).sin();
39            let sig2 = 0.7 * (2.0 * PI * 120.0 * t).sin();
40            let sig3 = 0.3 * (2.0 * PI * 300.0 * t).sin();
41
42            // Ruído simulado (pseudo-aleatório)
43            let noise = 0.1 * ((t * 12345.0).sin() + (t * 67890.0).cos());
44
45            sig1 + sig2 + sig3 + noise
46        })
47        .collect();
48
49    // Aplica janela de Blackman-Harris (alta atenuação)
50    println!("Aplicando janela de Blackman-Harris...");
51    let window_fn = window::blackman_harris::<f64>(n);
52    let windowed = window::apply(&signal, &window_fn);
53
54    // Converte para complexo
55    let complex_signal: Vec<Complex<f64>> = windowed
56        .iter()
57        .map(|&x| Complex::new(x, 0.0))
58        .collect();
59
60    // Cria planner e executa FFT
61    println!("Executando FFT...");
62    let planner = FftPlanner::new(n, false).unwrap();
63    let spectrum = planner.process(&complex_signal).unwrap();
64
65    // Analisa apenas metade do espectro (frequências positivas)
66    println!("\n=== Espectro de Potência ===\n");
67
68    let mut peaks = Vec::new();
69
70    for k in 0..n/2 {
71        let frequency = k as f64 * sample_rate / n as f64;
72        let magnitude = spectrum[k].norm();
73        let power_db = 20.0 * magnitude.log10();
74
75        // Detecta picos significativos (> -20 dB)
76        if power_db > -20.0 {
77            peaks.push((frequency, power_db));
78        }
79
80        // Imprime algumas frequências representativas
81        if k < 10 || (k % 50 == 0 && k < n/2) {
82            println!("  {:6.1} Hz: {:7.2} dB", frequency, power_db);
83        }
84    }
85
86    // Reporta picos detectados
87    println!("\n=== Picos Detectados (> -20 dB) ===\n");
88    peaks.sort_by(|a, b| b.1.partial_cmp(&a.1).unwrap());
89
90    for (i, (freq, db)) in peaks.iter().take(10).enumerate() {
91        println!("  {}. {:6.1} Hz: {:7.2} dB", i+1, freq, db);
92    }
93
94    // Teste de reversibilidade
95    println!("\n=== Teste de Reversibilidade ===\n");
96    let ifft_planner = FftPlanner::new(n, true).unwrap();
97    let recovered = ifft_planner.process(&spectrum).unwrap();
98
99    // Calcula erro RMS
100    let rms_error: f64 = complex_signal.iter()
101        .zip(recovered.iter())
102        .map(|(orig, rec)| (orig.re - rec.re).powi(2))
103        .sum::<f64>()
104        .sqrt() / (n as f64);
105
106    println!("Erro RMS (IFFT roundtrip): {:.2e}", rms_error);
107
108    if rms_error < 1e-10 {
109        println!("✓ Reversibilidade perfeita!");
110    }
111
112    // Demonstra diferentes janelas
113    println!("\n=== Comparação de Janelas ===\n");
114
115    let windows = vec![
116        ("Retangular", window::rectangular::<f64>(128)),
117        ("Hamming", window::hamming::<f64>(128)),
118        ("Hann", window::hann::<f64>(128)),
119        ("Blackman", window::blackman::<f64>(128)),
120    ];
121
122    for (name, win) in windows {
123        let sum: f64 = win.iter().sum();
124        let max = win.iter().cloned().fold(0.0, f64::max);
125        let min = win.iter().cloned().fold(1.0, f64::min);
126
127        println!("  {:<15} - soma: {:.2}, min: {:.3}, max: {:.3}",
128            name, sum, min, max);
129    }
130
131    println!("\n=== Análise Completa ===");
132}
spectral_analysis/spectral_analysis.rs

spectral_analysis/
spectral_analysis.rs
