//! Regime ADX Features
//!
//! This module implements ADX (Average Directional Index) and directional indicators
//! for regime detection. ADX measures trend strength using Wilder's smoothing method.
//!
//! **Wave D - Agent D14: ADX & Directional Indicators**
//! - Feature Indices: 211-215 (5 features)
//! - ADX (index 211): Trend strength (0-100)
//! - +DI (index 212): Positive directional indicator
//! - -DI (index 213): Negative directional indicator
//! - DX (index 214): Directional Movement Index (before ADX smoothing)
//! - ATR (index 215): Average True Range (volatility measure)
//!
//! ## Algorithm
//! 1. Calculate True Range (TR), +DM, -DM for each bar
//! 2. Apply Wilder's smoothing (similar to EMA with alpha=1/period)
//! 3. Compute +DI = 100 * smoothed_+DM / `smoothed_TR`
//! 4. Compute -DI = 100 * smoothed_-DM / `smoothed_TR`
//! 5. Compute DX = 100 * |+DI - -DI| / (+DI + -DI)
//! 6. Compute ADX = Wilder's smooth of DX
//!
//! ## Initialization Period
//! - First `period` bars: Simple moving average (SMA)
//! - After `period` bars: Wilder's smoothing (exponential-like)
//! - ADX initialization: Additional `period` bars after DX is ready
//! - Total warm-up: 2 * period bars (default: 28 bars)
//!
//! ## Performance Targets
//! - Per-feature calculation: <10μs
//! - Memory footprint: <200 bytes per symbol
//! - Cache-friendly: Sequential access patterns

pub use crate::OHLCVBar;

/// ADX Feature Extractor with Wilder's Smoothing
///
/// Maintains state for incremental ADX calculation using Wilder's smoothing method.
/// Requires 2 * period bars for full initialization (default: 28 bars).
#[derive(Debug)]
pub struct RegimeADXFeatures {
    /// Wilder's smoothing constant (1/period)
    alpha: f64,

    // Smoothed state variables
    /// Smoothed Average True Range
    atr: Option<f64>,
    /// Smoothed +DM (Positive Directional Movement)
    plus_dm_smooth: Option<f64>,
    /// Smoothed -DM (Negative Directional Movement)
    minus_dm_smooth: Option<f64>,
    /// Smoothed ADX value
    adx: Option<f64>,

    // Cache for feature output
    /// +DI value (0-100)
    plus_di: f64,
    /// -DI value (0-100)
    minus_di: f64,
    /// DX value (0-100)
    dx: f64,

    /// Previous bar for differential calculations
    prev_bar: Option<OHLCVBar>,
    /// Bar count for warmup tracking
    bar_count: usize,
}

impl RegimeADXFeatures {
    /// Create a new ADX feature extractor
    ///
    /// # Arguments
    /// * `period` - Smoothing period (default: 14, typical range: 7-28)
    ///
    /// # Returns
    /// A new `RegimeADXFeatures` instance with zero-initialized state
    ///
    /// # Example
    /// ```
    /// use ml::features::regime_adx::RegimeADXFeatures;
    ///
    /// let adx = RegimeADXFeatures::new(14);
    /// ```
    pub fn new(period: usize) -> Self {
        assert!(period > 0, "Period must be positive");

        Self {
            alpha: 1.0 / period as f64,
            atr: None,
            plus_dm_smooth: None,
            minus_dm_smooth: None,
            adx: None,
            plus_di: 0.0,
            minus_di: 0.0,
            dx: 0.0,
            prev_bar: None,
            bar_count: 0,
        }
    }

    /// Update ADX features with a new bar
    ///
    /// Returns 5 features:
    /// - `[0]`: ADX (0-100, trend strength)
    /// - `[1]`: +DI (0-100, positive directional indicator)
    /// - `[2]`: -DI (0-100, negative directional indicator)
    /// - `[3]`: DX (0-100, directional movement index)
    /// - `[4]`: ATR (>0, average true range)
    ///
    /// During initialization period (first bar), returns zeros.
    ///
    /// # Arguments
    /// * `bar` - OHLCV bar to process
    ///
    /// # Returns
    /// Array of 5 ADX-related features
    ///
    /// # Example
    /// ```ignore
    /// let features = adx.update(&bar);
    /// // features[0] is ADX, features[1] is +DI, etc.
    /// ```
    pub fn update(&mut self, bar: &OHLCVBar) -> [f64; 5] {
        self.bar_count += 1;

        // Need previous bar for calculations
        if self.prev_bar.is_none() {
            self.prev_bar = Some(*bar);
            return [0.0; 5];
        }

        let Some(prev) = self.prev_bar.as_ref() else {
            return [0.0; 5];
        };

        // 1. Calculate True Range (TR)
        let tr = self.calculate_true_range(bar, prev);

        // 2. Calculate Directional Movements (+DM, -DM)
        let (plus_dm, minus_dm) = self.calculate_directional_movements(bar, prev);

        // 3. Update smoothed values using Wilder's EMA
        self.update_smoothed_values(tr, plus_dm, minus_dm);

        // 4. Calculate Directional Indicators (+DI, -DI)
        self.calculate_directional_indicators();

        // 5. Calculate DX (Directional Index)
        self.calculate_dx();

        // 6. Update ADX (smoothed DX)
        self.update_adx();

        // Store current bar as previous
        self.prev_bar = Some(*bar);

        // Return feature vector: [ADX, +DI, -DI, DX, ATR]
        [
            self.adx.unwrap_or(0.0),
            self.plus_di,
            self.minus_di,
            self.dx,
            self.atr.unwrap_or(0.0),
        ]
    }

    /// Get current bar count (for warmup tracking)
    pub const fn bar_count(&self) -> usize {
        self.bar_count
    }

    /// Get Wilder's alpha constant
    pub const fn get_alpha(&self) -> f64 {
        self.alpha
    }

    // ========================================================================
    // Private Implementation Methods
    // ========================================================================

    /// Calculate True Range: max(H-L, |H-C_prev|, |L-C_prev|)
    fn calculate_true_range(&self, bar: &OHLCVBar, prev: &OHLCVBar) -> f64 {
        // WAVE 8 AGENT 32 FIX: Validate inputs are finite before arithmetic
        // If any input is NaN/Inf, return 0.0 to prevent NaN propagation
        if !bar.high.is_finite()
            || !bar.low.is_finite()
            || !bar.close.is_finite()
            || !prev.close.is_finite()
        {
            return 0.0;
        }

        let hl = bar.high - bar.low;
        let hc = (bar.high - prev.close).abs();
        let lc = (bar.low - prev.close).abs();
        let tr = hl.max(hc).max(lc);

        // Ensure TR is finite and non-negative
        if tr.is_finite() && tr >= 0.0 {
            tr
        } else {
            0.0
        }
    }

    /// Calculate +DM and -DM using Wilder's rules
    ///
    /// +DM = max(0, H - `H_prev`) if `high_diff` > `low_diff` and `high_diff` > 0
    /// -DM = max(0, `L_prev` - L) if `low_diff` > `high_diff` and `low_diff` > 0
    fn calculate_directional_movements(&self, bar: &OHLCVBar, prev: &OHLCVBar) -> (f64, f64) {
        // WAVE 8 AGENT 32 FIX: Validate inputs are finite before arithmetic
        // If any input is NaN/Inf, return (0.0, 0.0) to prevent NaN propagation
        if !bar.high.is_finite()
            || !bar.low.is_finite()
            || !prev.high.is_finite()
            || !prev.low.is_finite()
        {
            return (0.0, 0.0);
        }

        let high_diff = bar.high - prev.high;
        let low_diff = prev.low - bar.low;

        // Additional check: verify computed diffs are finite
        if !high_diff.is_finite() || !low_diff.is_finite() {
            return (0.0, 0.0);
        }

        let plus_dm = if high_diff > low_diff && high_diff > 0.0 {
            high_diff
        } else {
            0.0
        };

        let minus_dm = if low_diff > high_diff && low_diff > 0.0 {
            low_diff
        } else {
            0.0
        };

        (plus_dm, minus_dm)
    }

    /// Update smoothed values using Wilder's EMA: `S_new` = `S_old` × (1-α) + value × α
    fn update_smoothed_values(&mut self, tr: f64, plus_dm: f64, minus_dm: f64) {
        // ATR smoothing
        let new_atr = match self.atr {
            Some(prev) => prev * (1.0 - self.alpha) + tr * self.alpha,
            None => tr, // Initialize with first TR
        };
        self.atr = Some(if new_atr.is_finite() && new_atr >= 0.0 {
            new_atr
        } else {
            0.0
        });

        // +DM smoothing
        let new_plus_dm = match self.plus_dm_smooth {
            Some(prev) => prev * (1.0 - self.alpha) + plus_dm * self.alpha,
            None => plus_dm,
        };
        self.plus_dm_smooth = Some(if new_plus_dm.is_finite() && new_plus_dm >= 0.0 {
            new_plus_dm
        } else {
            0.0
        });

        // -DM smoothing
        let new_minus_dm = match self.minus_dm_smooth {
            Some(prev) => prev * (1.0 - self.alpha) + minus_dm * self.alpha,
            None => minus_dm,
        };
        self.minus_dm_smooth = Some(if new_minus_dm.is_finite() && new_minus_dm >= 0.0 {
            new_minus_dm
        } else {
            0.0
        });
    }

    /// Calculate +DI and -DI: DI = (`DM_smooth` / ATR) × 100
    fn calculate_directional_indicators(&mut self) {
        let atr_val = self.atr.unwrap_or(1.0);

        if atr_val > 1e-8 {
            self.plus_di = (self.plus_dm_smooth.unwrap_or(0.0) / atr_val) * 100.0;
            self.minus_di = (self.minus_dm_smooth.unwrap_or(0.0) / atr_val) * 100.0;
        } else {
            self.plus_di = 0.0;
            self.minus_di = 0.0;
        }

        // Ensure values are finite before clamping
        if !self.plus_di.is_finite() {
            self.plus_di = 0.0;
        }
        if !self.minus_di.is_finite() {
            self.minus_di = 0.0;
        }

        // Clamp to valid range [0, 100]
        self.plus_di = self.plus_di.clamp(0.0, 100.0);
        self.minus_di = self.minus_di.clamp(0.0, 100.0);
    }

    /// Calculate DX: |+DI - -DI| / (+DI + -DI) × 100
    fn calculate_dx(&mut self) {
        let di_sum = self.plus_di + self.minus_di;

        self.dx = if di_sum > 1e-8 {
            ((self.plus_di - self.minus_di).abs() / di_sum) * 100.0
        } else {
            0.0
        };

        // Ensure DX is finite
        if !self.dx.is_finite() {
            self.dx = 0.0;
        }

        // Clamp to valid range [0, 100]
        self.dx = self.dx.clamp(0.0, 100.0);
    }

    /// Update ADX using Wilder's smoothing of DX
    fn update_adx(&mut self) {
        self.adx = Some(match self.adx {
            Some(prev) => prev * (1.0 - self.alpha) + self.dx * self.alpha,
            None => self.dx, // Initialize with first DX
        });

        // Ensure ADX is finite and clamp to valid range [0, 100]
        if let Some(adx_val) = self.adx {
            if adx_val.is_finite() {
                self.adx = Some(adx_val.clamp(0.0, 100.0));
            } else {
                self.adx = Some(0.0);
            }
        }
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    /// Helper to create a test bar
    fn create_test_bar(open: f64, high: f64, low: f64, close: f64, volume: f64) -> OHLCVBar {
        OHLCVBar {
            timestamp: chrono::DateTime::<chrono::Utc>::from_timestamp(0, 0).unwrap(),
            open,
            high,
            low,
            close,
            volume,
        }
    }

    #[test]
    fn test_new_initialization() {
        let adx = RegimeADXFeatures::new(14);
        assert!((adx.alpha - 1.0 / 14.0).abs() < 1e-10);
        assert!(adx.atr.is_none());
        assert!(adx.adx.is_none());
        assert!(adx.prev_bar.is_none());
        assert_eq!(adx.bar_count, 0);
    }

    #[test]
    fn test_update_returns_zeros_initially() {
        let mut adx = RegimeADXFeatures::new(14);
        let bar = create_test_bar(100.0, 102.0, 99.0, 101.0, 1000.0);

        let features = adx.update(&bar);
        assert_eq!(features, [0.0; 5]);
    }

    #[test]
    fn test_custom_period() {
        let adx = RegimeADXFeatures::new(20);
        assert!((adx.alpha - 1.0 / 20.0).abs() < 1e-10);
    }

    /// WAVE 8 AGENT 32: Test ADX handles NaN inputs gracefully (no NaN propagation)
    #[test]
    fn test_adx_handles_nan_inputs() {
        let mut adx = RegimeADXFeatures::new(14);

        // First valid bar
        let bar1 = create_test_bar(100.0, 102.0, 98.0, 101.0, 1000.0);
        let features1 = adx.update(&bar1);
        assert_eq!(features1, [0.0; 5]); // First bar returns zeros

        // Second bar with NaN high (simulates corrupted DBN data)
        let bar2_nan = OHLCVBar {
            timestamp: chrono::DateTime::<chrono::Utc>::from_timestamp(0, 0).unwrap(),
            open: 101.0,
            high: f64::NAN, // NaN input - triggers the fix
            low: 99.0,
            close: 100.0,
            volume: 1000.0,
        };
        let features2 = adx.update(&bar2_nan);

        // Should handle NaN gracefully - return finite values (zeros or valid)
        assert!(
            features2[0].is_finite(),
            "ADX should be finite, got: {}",
            features2[0]
        );
        assert!(
            features2[1].is_finite(),
            "+DI should be finite, got: {}",
            features2[1]
        );
        assert!(
            features2[2].is_finite(),
            "-DI should be finite, got: {}",
            features2[2]
        );
        assert!(
            features2[3].is_finite(),
            "DX should be finite, got: {}",
            features2[3]
        );
        assert!(
            features2[4].is_finite(),
            "ATR should be finite, got: {}",
            features2[4]
        );
    }

    /// WAVE 8 AGENT 32: Test ADX handles Inf inputs gracefully
    #[test]
    fn test_adx_handles_inf_inputs() {
        let mut adx = RegimeADXFeatures::new(14);

        // First valid bar
        let bar1 = create_test_bar(100.0, 102.0, 98.0, 101.0, 1000.0);
        adx.update(&bar1);

        // Second bar with Inf low (simulates price anomaly)
        let bar2_inf = OHLCVBar {
            timestamp: chrono::DateTime::<chrono::Utc>::from_timestamp(0, 0).unwrap(),
            open: 101.0,
            high: 103.0,
            low: f64::INFINITY, // Inf input - triggers the fix
            close: 100.0,
            volume: 1000.0,
        };
        let features2 = adx.update(&bar2_inf);

        // Should handle Inf gracefully
        assert!(features2[0].is_finite(), "ADX should be finite");
        assert!(features2[1].is_finite(), "+DI should be finite");
        assert!(features2[2].is_finite(), "-DI should be finite");
        assert!(features2[3].is_finite(), "DX should be finite");
        assert!(features2[4].is_finite(), "ATR should be finite");
    }

    /// WAVE 8 AGENT 32: Test ADX with multiple consecutive NaN bars
    #[test]
    fn test_adx_multiple_nan_bars() {
        let mut adx = RegimeADXFeatures::new(14);

        // Feed 10 bars with various NaN values
        let epoch = chrono::DateTime::<chrono::Utc>::from_timestamp(0, 0).unwrap();
        for i in 0..10 {
            let bar = OHLCVBar {
                timestamp: epoch,
                open: if i % 2 == 0 { 100.0 } else { f64::NAN },
                high: if i % 3 == 0 { 102.0 } else { f64::NAN },
                low: if i % 4 == 0 { 98.0 } else { f64::NAN },
                close: if i % 5 == 0 { 101.0 } else { f64::NAN },
                volume: 1000.0,
            };
            let features = adx.update(&bar);

            // All features should remain finite
            for (idx, &feat) in features.iter().enumerate() {
                assert!(
                    feat.is_finite(),
                    "Feature {} should be finite at bar {}, got: {}",
                    idx,
                    i,
                    feat
                );
            }
        }
    }
}
