oxigdal_vrt/
band.rs

1//! VRT virtual band configuration
2
3use crate::error::{Result, VrtError};
4use crate::source::VrtSource;
5use oxigdal_core::types::{ColorInterpretation, NoDataValue, RasterDataType};
6use serde::{Deserialize, Serialize};
7
8/// VRT band configuration
9#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
10pub struct VrtBand {
11    /// Band number (1-based)
12    pub band: usize,
13    /// Data type
14    pub data_type: RasterDataType,
15    /// Color interpretation
16    pub color_interp: ColorInterpretation,
17    /// NoData value
18    pub nodata: NoDataValue,
19    /// Sources for this band
20    pub sources: Vec<VrtSource>,
21    /// Block size (tile dimensions)
22    pub block_size: Option<(u32, u32)>,
23    /// Pixel function (for on-the-fly computation)
24    pub pixel_function: Option<PixelFunction>,
25    /// Offset for scaling
26    pub offset: Option<f64>,
27    /// Scale factor
28    pub scale: Option<f64>,
29    /// Color table
30    pub color_table: Option<ColorTable>,
31}
32
33impl VrtBand {
34    /// Creates a new VRT band
35    pub fn new(band: usize, data_type: RasterDataType) -> Self {
36        Self {
37            band,
38            data_type,
39            color_interp: ColorInterpretation::Undefined,
40            nodata: NoDataValue::None,
41            sources: Vec::new(),
42            block_size: None,
43            pixel_function: None,
44            offset: None,
45            scale: None,
46            color_table: None,
47        }
48    }
49
50    /// Creates a simple band with a single source
51    pub fn simple(band: usize, data_type: RasterDataType, source: VrtSource) -> Self {
52        Self {
53            band,
54            data_type,
55            color_interp: ColorInterpretation::Undefined,
56            nodata: source.nodata.unwrap_or(NoDataValue::None),
57            sources: vec![source],
58            block_size: None,
59            pixel_function: None,
60            offset: None,
61            scale: None,
62            color_table: None,
63        }
64    }
65
66    /// Adds a source to this band
67    pub fn add_source(&mut self, source: VrtSource) {
68        self.sources.push(source);
69    }
70
71    /// Sets the color interpretation
72    pub fn with_color_interp(mut self, color_interp: ColorInterpretation) -> Self {
73        self.color_interp = color_interp;
74        self
75    }
76
77    /// Sets the NoData value
78    pub fn with_nodata(mut self, nodata: NoDataValue) -> Self {
79        self.nodata = nodata;
80        self
81    }
82
83    /// Sets the block size
84    pub fn with_block_size(mut self, width: u32, height: u32) -> Self {
85        self.block_size = Some((width, height));
86        self
87    }
88
89    /// Sets the pixel function
90    pub fn with_pixel_function(mut self, function: PixelFunction) -> Self {
91        self.pixel_function = Some(function);
92        self
93    }
94
95    /// Sets the offset and scale
96    pub fn with_scaling(mut self, offset: f64, scale: f64) -> Self {
97        self.offset = Some(offset);
98        self.scale = Some(scale);
99        self
100    }
101
102    /// Sets the color table
103    pub fn with_color_table(mut self, color_table: ColorTable) -> Self {
104        self.color_table = Some(color_table);
105        self
106    }
107
108    /// Validates the band configuration
109    ///
110    /// # Errors
111    /// Returns an error if the band is invalid
112    pub fn validate(&self) -> Result<()> {
113        if self.band == 0 {
114            return Err(VrtError::invalid_band("Band number must be >= 1"));
115        }
116
117        if self.sources.is_empty() && self.pixel_function.is_none() {
118            return Err(VrtError::invalid_band(
119                "Band must have at least one source or a pixel function",
120            ));
121        }
122
123        // Validate all sources
124        for source in &self.sources {
125            source.validate()?;
126        }
127
128        // Validate pixel function if present
129        if let Some(ref func) = self.pixel_function {
130            func.validate(&self.sources)?;
131        }
132
133        Ok(())
134    }
135
136    /// Checks if this band has multiple sources
137    pub fn has_multiple_sources(&self) -> bool {
138        self.sources.len() > 1
139    }
140
141    /// Checks if this band uses a pixel function
142    pub fn uses_pixel_function(&self) -> bool {
143        self.pixel_function.is_some()
144    }
145
146    /// Applies scaling to a value
147    pub fn apply_scaling(&self, value: f64) -> f64 {
148        let scaled = if let Some(scale) = self.scale {
149            value * scale
150        } else {
151            value
152        };
153
154        if let Some(offset) = self.offset {
155            scaled + offset
156        } else {
157            scaled
158        }
159    }
160}
161
162/// Pixel function for on-the-fly computation
163#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
164pub enum PixelFunction {
165    /// Average of all sources
166    Average,
167    /// Minimum of all sources
168    Min,
169    /// Maximum of all sources
170    Max,
171    /// Sum of all sources
172    Sum,
173    /// First valid (non-NoData) value
174    FirstValid,
175    /// Last valid (non-NoData) value
176    LastValid,
177    /// Weighted average (requires weights)
178    WeightedAverage {
179        /// Weights for each source
180        weights: Vec<f64>,
181    },
182    /// NDVI: (NIR - Red) / (NIR + Red)
183    /// Requires exactly 2 sources: [Red, NIR]
184    Ndvi,
185    /// Enhanced Vegetation Index: 2.5 * (NIR - Red) / (NIR + 6*Red - 7.5*Blue + 1)
186    /// Requires exactly 3 sources: [Red, NIR, Blue]
187    Evi,
188    /// Normalized Difference Water Index: (Green - NIR) / (Green + NIR)
189    /// Requires exactly 2 sources: [Green, NIR]
190    Ndwi,
191    /// Band math expression
192    /// Supports operations: +, -, *, /, sqrt, pow, abs, min, max
193    /// Variables are named as B1, B2, B3, etc. corresponding to source bands
194    BandMath {
195        /// Expression string (e.g., "(B1 + B2) / 2", "sqrt(B1 * B2)")
196        expression: String,
197    },
198    /// Lookup table transformation
199    /// Maps input values to output values
200    LookupTable {
201        /// Lookup table: vec of (input_value, output_value) pairs
202        table: Vec<(f64, f64)>,
203        /// Interpolation method: "nearest", "linear"
204        interpolation: String,
205    },
206    /// Conditional logic: if condition then value_if_true else value_if_false
207    /// Condition format: "B1 > 0.5", "B1 >= B2", etc.
208    Conditional {
209        /// Condition expression
210        condition: String,
211        /// Value if condition is true (can be a constant or expression)
212        value_if_true: String,
213        /// Value if condition is false (can be a constant or expression)
214        value_if_false: String,
215    },
216    /// Multiply source values
217    Multiply,
218    /// Divide first source by second (handles division by zero)
219    Divide,
220    /// Square root of source value
221    SquareRoot,
222    /// Absolute value of source value
223    Absolute,
224    /// Custom function (not yet implemented)
225    Custom {
226        /// Function name
227        name: String,
228    },
229}
230
231impl PixelFunction {
232    /// Validates the pixel function against sources
233    ///
234    /// # Errors
235    /// Returns an error if the function is invalid for the given sources
236    pub fn validate(&self, sources: &[VrtSource]) -> Result<()> {
237        match self {
238            Self::WeightedAverage { weights } => {
239                if weights.len() != sources.len() {
240                    return Err(VrtError::invalid_band(format!(
241                        "WeightedAverage requires {} weights, got {}",
242                        sources.len(),
243                        weights.len()
244                    )));
245                }
246
247                // Check that weights sum to approximately 1.0
248                let sum: f64 = weights.iter().sum();
249                if (sum - 1.0).abs() > 0.001 {
250                    return Err(VrtError::invalid_band(format!(
251                        "Weights should sum to 1.0, got {}",
252                        sum
253                    )));
254                }
255            }
256            Self::Ndvi | Self::Ndwi if sources.len() != 2 => {
257                return Err(VrtError::invalid_band(format!(
258                    "{:?} requires exactly 2 sources, got {}",
259                    self,
260                    sources.len()
261                )));
262            }
263            Self::Evi if sources.len() != 3 => {
264                return Err(VrtError::invalid_band(format!(
265                    "EVI requires exactly 3 sources, got {}",
266                    sources.len()
267                )));
268            }
269            Self::BandMath { expression } if expression.trim().is_empty() => {
270                return Err(VrtError::invalid_band(
271                    "BandMath expression cannot be empty",
272                ));
273            }
274            Self::LookupTable { table, .. } if table.is_empty() => {
275                return Err(VrtError::invalid_band("LookupTable cannot be empty"));
276            }
277            Self::Conditional {
278                condition,
279                value_if_true,
280                value_if_false,
281            } => {
282                if condition.trim().is_empty() {
283                    return Err(VrtError::invalid_band(
284                        "Conditional condition cannot be empty",
285                    ));
286                }
287                if value_if_true.trim().is_empty() || value_if_false.trim().is_empty() {
288                    return Err(VrtError::invalid_band("Conditional values cannot be empty"));
289                }
290            }
291            Self::Divide | Self::Multiply if sources.len() < 2 => {
292                return Err(VrtError::invalid_band(format!(
293                    "{:?} requires at least 2 sources, got {}",
294                    self,
295                    sources.len()
296                )));
297            }
298            Self::SquareRoot | Self::Absolute if sources.is_empty() => {
299                return Err(VrtError::invalid_band(format!(
300                    "{:?} requires at least 1 source",
301                    self
302                )));
303            }
304            Self::Custom { name } => {
305                return Err(VrtError::InvalidPixelFunction {
306                    function: name.clone(),
307                });
308            }
309            _ => {}
310        }
311        Ok(())
312    }
313
314    /// Applies the pixel function to a set of values
315    ///
316    /// # Errors
317    /// Returns an error if the function cannot be applied
318    pub fn apply(&self, values: &[Option<f64>]) -> Result<Option<f64>> {
319        match self {
320            Self::Average => {
321                let valid: Vec<f64> = values.iter().filter_map(|v| *v).collect();
322                if valid.is_empty() {
323                    Ok(None)
324                } else {
325                    Ok(Some(valid.iter().sum::<f64>() / valid.len() as f64))
326                }
327            }
328            Self::Min => Ok(values
329                .iter()
330                .filter_map(|v| *v)
331                .min_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal))),
332            Self::Max => Ok(values
333                .iter()
334                .filter_map(|v| *v)
335                .max_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal))),
336            Self::Sum => {
337                let valid: Vec<f64> = values.iter().filter_map(|v| *v).collect();
338                if valid.is_empty() {
339                    Ok(None)
340                } else {
341                    Ok(Some(valid.iter().sum()))
342                }
343            }
344            Self::FirstValid => Ok(values.iter().find_map(|v| *v)),
345            Self::LastValid => Ok(values.iter().rev().find_map(|v| *v)),
346            Self::WeightedAverage { weights } => {
347                if weights.len() != values.len() {
348                    return Err(VrtError::invalid_band("Weight count mismatch"));
349                }
350
351                let mut sum = 0.0;
352                let mut weight_sum = 0.0;
353
354                for (value, weight) in values.iter().zip(weights.iter()) {
355                    if let Some(v) = value {
356                        sum += v * weight;
357                        weight_sum += weight;
358                    }
359                }
360
361                if weight_sum > 0.0 {
362                    Ok(Some(sum / weight_sum))
363                } else {
364                    Ok(None)
365                }
366            }
367            Self::Ndvi => {
368                // NDVI = (NIR - Red) / (NIR + Red)
369                if values.len() != 2 {
370                    return Err(VrtError::invalid_band("NDVI requires exactly 2 values"));
371                }
372                match (values[0], values[1]) {
373                    (Some(red), Some(nir)) => {
374                        let denominator = nir + red;
375                        if denominator.abs() < f64::EPSILON {
376                            Ok(None) // Avoid division by zero
377                        } else {
378                            Ok(Some((nir - red) / denominator))
379                        }
380                    }
381                    _ => Ok(None),
382                }
383            }
384            Self::Evi => {
385                // EVI = 2.5 * (NIR - Red) / (NIR + 6*Red - 7.5*Blue + 1)
386                if values.len() != 3 {
387                    return Err(VrtError::invalid_band("EVI requires exactly 3 values"));
388                }
389                match (values[0], values[1], values[2]) {
390                    (Some(red), Some(nir), Some(blue)) => {
391                        let denominator = nir + 6.0 * red - 7.5 * blue + 1.0;
392                        if denominator.abs() < f64::EPSILON {
393                            Ok(None)
394                        } else {
395                            Ok(Some(2.5 * (nir - red) / denominator))
396                        }
397                    }
398                    _ => Ok(None),
399                }
400            }
401            Self::Ndwi => {
402                // NDWI = (Green - NIR) / (Green + NIR)
403                if values.len() != 2 {
404                    return Err(VrtError::invalid_band("NDWI requires exactly 2 values"));
405                }
406                match (values[0], values[1]) {
407                    (Some(green), Some(nir)) => {
408                        let denominator = green + nir;
409                        if denominator.abs() < f64::EPSILON {
410                            Ok(None)
411                        } else {
412                            Ok(Some((green - nir) / denominator))
413                        }
414                    }
415                    _ => Ok(None),
416                }
417            }
418            Self::BandMath { expression } => Self::evaluate_expression(expression, values),
419            Self::LookupTable {
420                table,
421                interpolation,
422            } => {
423                if values.is_empty() {
424                    return Ok(None);
425                }
426                if let Some(value) = values[0] {
427                    Self::apply_lookup_table(value, table, interpolation)
428                } else {
429                    Ok(None)
430                }
431            }
432            Self::Conditional {
433                condition,
434                value_if_true,
435                value_if_false,
436            } => Self::evaluate_conditional(condition, value_if_true, value_if_false, values),
437            Self::Multiply => {
438                let valid: Vec<f64> = values.iter().filter_map(|v| *v).collect();
439                if valid.is_empty() {
440                    Ok(None)
441                } else {
442                    Ok(Some(valid.iter().product()))
443                }
444            }
445            Self::Divide => {
446                if values.len() < 2 {
447                    return Err(VrtError::invalid_band("Divide requires at least 2 values"));
448                }
449                match (values[0], values[1]) {
450                    (Some(numerator), Some(denominator)) => {
451                        if denominator.abs() < f64::EPSILON {
452                            Ok(None) // Avoid division by zero
453                        } else {
454                            Ok(Some(numerator / denominator))
455                        }
456                    }
457                    _ => Ok(None),
458                }
459            }
460            Self::SquareRoot => {
461                if values.is_empty() {
462                    return Ok(None);
463                }
464                values[0].map_or(Ok(None), |v| {
465                    if v < 0.0 {
466                        Ok(None) // Negative values have no real square root
467                    } else {
468                        Ok(Some(v.sqrt()))
469                    }
470                })
471            }
472            Self::Absolute => {
473                if values.is_empty() {
474                    return Ok(None);
475                }
476                Ok(values[0].map(|v| v.abs()))
477            }
478            Self::Custom { name } => Err(VrtError::InvalidPixelFunction {
479                function: name.clone(),
480            }),
481        }
482    }
483
484    /// Evaluates a band math expression
485    fn evaluate_expression(expression: &str, values: &[Option<f64>]) -> Result<Option<f64>> {
486        // Simple expression evaluator for basic band math
487        // Supports: +, -, *, /, sqrt, pow, abs, min, max
488        // Variables: B1, B2, B3, etc.
489
490        // If any band is NoData, the result is NoData (preserves prior semantics).
491        if values.iter().any(Option::is_none) {
492            return Ok(None);
493        }
494
495        // Substitute band variables (B1, B2, ... B10, B11, ...) with their
496        // values. A naive `String::replace("B1", ...)` would corrupt "B10"/"B11"
497        // by matching the "B1" prefix, so we scan for a `B` followed by a full
498        // run of digits and substitute the *complete* index token atomically.
499        let expr = Self::substitute_band_variables(expression, values);
500
501        // Basic expression evaluation (simplified)
502        // For production, consider using a proper expression parser like `evalexpr`
503        match Self::simple_eval(&expr) {
504            Ok(result) => Ok(Some(result)),
505            Err(_) => Err(VrtError::invalid_band(format!(
506                "Failed to evaluate expression: {}",
507                expression
508            ))),
509        }
510    }
511
512    /// Substitutes band variables (`B1`, `B2`, ..., `B10`, `B11`, ...) in
513    /// `expression` with their numeric values from `values` (1-indexed).
514    ///
515    /// Unlike a plain `String::replace`, this matches each `B<n>` token as a
516    /// complete unit — a `B` followed by a maximal run of ASCII digits — so
517    /// `B1` never matches inside `B10`/`B11`. Tokens whose index is out of range
518    /// (no corresponding band) are left untouched. Callers are expected to have
519    /// already rejected NoData bands; a `None` value is treated as `0.0` here
520    /// only as a defensive fallback.
521    fn substitute_band_variables(expression: &str, values: &[Option<f64>]) -> String {
522        let chars: Vec<char> = expression.chars().collect();
523        let mut out = String::with_capacity(expression.len());
524        let mut i = 0;
525        while i < chars.len() {
526            if chars[i] == 'B' && i + 1 < chars.len() && chars[i + 1].is_ascii_digit() {
527                // Consume the full run of digits into a band index.
528                let mut j = i + 1;
529                let mut index: usize = 0;
530                while j < chars.len() && chars[j].is_ascii_digit() {
531                    index = index
532                        .saturating_mul(10)
533                        .saturating_add((chars[j] as u8 - b'0') as usize);
534                    j += 1;
535                }
536                if index >= 1 && index <= values.len() {
537                    let value = values[index - 1].unwrap_or(0.0);
538                    out.push_str(&value.to_string());
539                } else {
540                    // Not a valid band reference: keep the original token.
541                    out.extend(&chars[i..j]);
542                }
543                i = j;
544            } else {
545                out.push(chars[i]);
546                i += 1;
547            }
548        }
549        out
550    }
551
552    /// Simple expression evaluator (basic implementation)
553    fn simple_eval(expr: &str) -> Result<f64> {
554        let expr = expr.trim();
555
556        // Try to parse as number first
557        if let Ok(num) = expr.parse::<f64>() {
558            return Ok(num);
559        }
560
561        // Handle sqrt
562        if expr.starts_with("sqrt(") && expr.ends_with(')') {
563            let inner = &expr[5..expr.len() - 1];
564            let val = Self::simple_eval(inner)?;
565            if val < 0.0 {
566                return Err(VrtError::invalid_band("Square root of negative number"));
567            }
568            return Ok(val.sqrt());
569        }
570
571        // Handle abs
572        if expr.starts_with("abs(") && expr.ends_with(')') {
573            let inner = &expr[4..expr.len() - 1];
574            let val = Self::simple_eval(inner)?;
575            return Ok(val.abs());
576        }
577
578        // Handle parentheses
579        if expr.starts_with('(') && expr.ends_with(')') {
580            // Check if these are balanced outer parens
581            let inner = &expr[1..expr.len() - 1];
582            let mut depth = 0;
583            let mut is_outer = true;
584            for ch in inner.chars() {
585                if ch == '(' {
586                    depth += 1;
587                } else if ch == ')' {
588                    depth -= 1;
589                    if depth < 0 {
590                        is_outer = false;
591                        break;
592                    }
593                }
594            }
595            if is_outer && depth == 0 {
596                return Self::simple_eval(inner);
597            }
598        }
599
600        // Handle binary operations (search from right to left for left-to-right evaluation)
601        // Process + and - first (lower precedence), then * and / (higher precedence)
602        for op in &['+', '-'] {
603            let mut depth = 0;
604            for (i, ch) in expr.char_indices().rev() {
605                if ch == ')' {
606                    depth += 1;
607                } else if ch == '(' {
608                    depth -= 1;
609                } else if depth == 0 && ch == *op && i > 0 && i < expr.len() - 1 {
610                    let left = Self::simple_eval(&expr[..i])?;
611                    let right = Self::simple_eval(&expr[i + 1..])?;
612                    return match op {
613                        '+' => Ok(left + right),
614                        '-' => Ok(left - right),
615                        _ => unreachable!(),
616                    };
617                }
618            }
619        }
620
621        for op in &['*', '/'] {
622            let mut depth = 0;
623            for (i, ch) in expr.char_indices().rev() {
624                if ch == ')' {
625                    depth += 1;
626                } else if ch == '(' {
627                    depth -= 1;
628                } else if depth == 0 && ch == *op && i > 0 && i < expr.len() - 1 {
629                    let left = Self::simple_eval(&expr[..i])?;
630                    let right = Self::simple_eval(&expr[i + 1..])?;
631                    return match op {
632                        '*' => Ok(left * right),
633                        '/' => {
634                            if right.abs() < f64::EPSILON {
635                                Err(VrtError::invalid_band("Division by zero"))
636                            } else {
637                                Ok(left / right)
638                            }
639                        }
640                        _ => unreachable!(),
641                    };
642                }
643            }
644        }
645
646        Err(VrtError::invalid_band(format!(
647            "Cannot parse expression: {}",
648            expr
649        )))
650    }
651
652    /// Applies lookup table transformation
653    fn apply_lookup_table(
654        value: f64,
655        table: &[(f64, f64)],
656        interpolation: &str,
657    ) -> Result<Option<f64>> {
658        if table.is_empty() {
659            return Ok(None);
660        }
661
662        match interpolation {
663            "nearest" => {
664                // Find nearest entry
665                let mut best_idx = 0;
666                let mut best_dist = (table[0].0 - value).abs();
667
668                for (i, (input, _)) in table.iter().enumerate() {
669                    let dist = (input - value).abs();
670                    if dist < best_dist {
671                        best_dist = dist;
672                        best_idx = i;
673                    }
674                }
675
676                Ok(Some(table[best_idx].1))
677            }
678            "linear" => {
679                // Linear interpolation
680                if value <= table[0].0 {
681                    return Ok(Some(table[0].1));
682                }
683                if value >= table[table.len() - 1].0 {
684                    return Ok(Some(table[table.len() - 1].1));
685                }
686
687                // Find surrounding points
688                for i in 0..table.len() - 1 {
689                    if value >= table[i].0 && value <= table[i + 1].0 {
690                        let x0 = table[i].0;
691                        let y0 = table[i].1;
692                        let x1 = table[i + 1].0;
693                        let y1 = table[i + 1].1;
694
695                        let t = (value - x0) / (x1 - x0);
696                        return Ok(Some(y0 + t * (y1 - y0)));
697                    }
698                }
699
700                Ok(Some(table[0].1))
701            }
702            _ => Err(VrtError::invalid_band(format!(
703                "Unknown interpolation method: {}",
704                interpolation
705            ))),
706        }
707    }
708
709    /// Evaluates conditional expression
710    fn evaluate_conditional(
711        condition: &str,
712        value_if_true: &str,
713        value_if_false: &str,
714        values: &[Option<f64>],
715    ) -> Result<Option<f64>> {
716        // Simple condition evaluator
717        // Supports: >, <, >=, <=, ==, !=
718
719        let cond_result = Self::evaluate_condition(condition, values)?;
720
721        let target_expr = if cond_result {
722            value_if_true
723        } else {
724            value_if_false
725        };
726
727        // Evaluate the target expression
728        Self::evaluate_expression(target_expr, values)
729    }
730
731    /// Evaluates a boolean condition
732    fn evaluate_condition(condition: &str, values: &[Option<f64>]) -> Result<bool> {
733        let condition = condition.trim();
734
735        // Try different comparison operators in order (longer ones first to avoid false matches)
736        let operators = [">=", "<=", "==", "!=", ">", "<"];
737
738        for op_str in &operators {
739            if let Some(pos) = condition.find(op_str) {
740                let left_expr = condition[..pos].trim();
741                let right_expr = condition[pos + op_str.len()..].trim();
742
743                let left_val = Self::evaluate_expression(left_expr, values)?
744                    .ok_or_else(|| VrtError::invalid_band("Left side of condition is NoData"))?;
745
746                let right_val = Self::evaluate_expression(right_expr, values)?
747                    .ok_or_else(|| VrtError::invalid_band("Right side of condition is NoData"))?;
748
749                let result = match *op_str {
750                    ">=" => left_val >= right_val,
751                    "<=" => left_val <= right_val,
752                    ">" => left_val > right_val,
753                    "<" => left_val < right_val,
754                    "==" => (left_val - right_val).abs() < f64::EPSILON,
755                    "!=" => (left_val - right_val).abs() >= f64::EPSILON,
756                    _ => {
757                        return Err(VrtError::invalid_band(format!(
758                            "Unknown operator: {}",
759                            op_str
760                        )));
761                    }
762                };
763
764                return Ok(result);
765            }
766        }
767
768        Err(VrtError::invalid_band(format!(
769            "Cannot parse condition: {}",
770            condition
771        )))
772    }
773}
774
775/// Color table entry
776#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
777pub struct ColorEntry {
778    /// Color value (index)
779    pub value: u16,
780    /// Red component (0-255)
781    pub r: u8,
782    /// Green component (0-255)
783    pub g: u8,
784    /// Blue component (0-255)
785    pub b: u8,
786    /// Alpha component (0-255)
787    pub a: u8,
788}
789
790impl ColorEntry {
791    /// Creates a new color entry
792    pub const fn new(value: u16, r: u8, g: u8, b: u8, a: u8) -> Self {
793        Self { value, r, g, b, a }
794    }
795
796    /// Creates an opaque color entry
797    pub const fn rgb(value: u16, r: u8, g: u8, b: u8) -> Self {
798        Self::new(value, r, g, b, 255)
799    }
800}
801
802/// Color table (palette)
803#[derive(Debug, Clone, PartialEq, Eq, Serialize, Deserialize)]
804pub struct ColorTable {
805    /// Color entries
806    pub entries: Vec<ColorEntry>,
807}
808
809impl ColorTable {
810    /// Creates a new empty color table
811    pub fn new() -> Self {
812        Self {
813            entries: Vec::new(),
814        }
815    }
816
817    /// Creates a color table with entries
818    pub fn with_entries(entries: Vec<ColorEntry>) -> Self {
819        Self { entries }
820    }
821
822    /// Adds a color entry
823    pub fn add_entry(&mut self, entry: ColorEntry) {
824        self.entries.push(entry);
825    }
826
827    /// Gets a color entry by value
828    pub fn get(&self, value: u16) -> Option<&ColorEntry> {
829        self.entries.iter().find(|e| e.value == value)
830    }
831}
832
833impl Default for ColorTable {
834    fn default() -> Self {
835        Self::new()
836    }
837}
838
839#[cfg(test)]
840mod tests {
841    use super::*;
842    use crate::source::SourceFilename;
843
844    #[test]
845    fn test_vrt_band_creation() {
846        let band = VrtBand::new(1, RasterDataType::UInt8);
847        assert_eq!(band.band, 1);
848        assert_eq!(band.data_type, RasterDataType::UInt8);
849    }
850
851    #[test]
852    fn test_vrt_band_validation() {
853        let source = VrtSource::new(SourceFilename::absolute("/test.tif"), 1);
854        let band = VrtBand::simple(1, RasterDataType::UInt8, source);
855        assert!(band.validate().is_ok());
856
857        let invalid_band = VrtBand::new(0, RasterDataType::UInt8);
858        assert!(invalid_band.validate().is_err());
859
860        let no_source_band = VrtBand::new(1, RasterDataType::UInt8);
861        assert!(no_source_band.validate().is_err());
862    }
863
864    #[test]
865    fn test_pixel_function_average() {
866        let func = PixelFunction::Average;
867        let values = vec![Some(1.0), Some(2.0), Some(3.0)];
868        let result = func.apply(&values);
869        assert!(result.is_ok());
870        assert_eq!(result.ok().flatten(), Some(2.0));
871
872        let with_none = vec![Some(1.0), None, Some(3.0)];
873        let result = func.apply(&with_none);
874        assert!(result.is_ok());
875        assert_eq!(result.ok().flatten(), Some(2.0));
876    }
877
878    #[test]
879    fn test_pixel_function_first_valid() {
880        let func = PixelFunction::FirstValid;
881        let values = vec![None, Some(2.0), Some(3.0)];
882        let result = func.apply(&values);
883        assert!(result.is_ok());
884        assert_eq!(result.ok().flatten(), Some(2.0));
885    }
886
887    #[test]
888    fn test_pixel_function_weighted_average() {
889        let func = PixelFunction::WeightedAverage {
890            weights: vec![0.5, 0.3, 0.2],
891        };
892        let values = vec![Some(10.0), Some(20.0), Some(30.0)];
893        let result = func.apply(&values);
894        assert!(result.is_ok());
895        // 10*0.5 + 20*0.3 + 30*0.2 = 5 + 6 + 6 = 17
896        assert_eq!(result.ok().flatten(), Some(17.0));
897    }
898
899    #[test]
900    fn test_band_scaling() {
901        let band = VrtBand::new(1, RasterDataType::UInt8).with_scaling(10.0, 2.0);
902
903        assert_eq!(band.apply_scaling(5.0), 20.0); // 5 * 2 + 10 = 20
904    }
905
906    #[test]
907    fn test_color_table() {
908        let mut table = ColorTable::new();
909        table.add_entry(ColorEntry::rgb(0, 255, 0, 0));
910        table.add_entry(ColorEntry::rgb(1, 0, 255, 0));
911        table.add_entry(ColorEntry::rgb(2, 0, 0, 255));
912
913        assert_eq!(table.entries.len(), 3);
914        assert_eq!(table.get(1).map(|e| e.g), Some(255));
915    }
916
917    /// Regression: BandMath expressions referencing 10+ bands must substitute
918    /// each `B<n>` token independently. The old `String::replace("B1", ...)`
919    /// corrupted `B10`/`B11` by matching the `B1` prefix.
920    #[test]
921    fn test_band_math_double_digit_bands() {
922        let func = PixelFunction::BandMath {
923            expression: "B1+B10+B11".to_string(),
924        };
925        // 11 bands: B1=1, B2..B9=0, B10=10, B11=11.
926        let mut values = vec![Some(0.0); 11];
927        values[0] = Some(1.0); // B1
928        values[9] = Some(10.0); // B10
929        values[10] = Some(11.0); // B11
930
931        let result = func.apply(&values).expect("band math should evaluate");
932        assert_eq!(result, Some(22.0), "expected 1 + 10 + 11 = 22");
933    }
934
935    /// A lone `B10` reference must resolve to band 10's value, never band 1's
936    /// value spliced into the token.
937    #[test]
938    fn test_band_math_b10_not_contaminated_by_b1() {
939        let func = PixelFunction::BandMath {
940            expression: "B10".to_string(),
941        };
942        let mut values = vec![Some(0.0); 10];
943        values[0] = Some(7.0); // B1 — distinctive, must NOT leak into B10
944        values[9] = Some(3.0); // B10
945
946        let result = func.apply(&values).expect("band math should evaluate");
947        assert_eq!(result, Some(3.0));
948    }
949
950    /// The Conditional pixel function evaluates its sub-expressions through the
951    /// same substitution path, so it must also handle 10+ band references.
952    #[test]
953    fn test_conditional_double_digit_bands() {
954        let func = PixelFunction::Conditional {
955            condition: "B11>B1".to_string(),
956            value_if_true: "B10".to_string(),
957            value_if_false: "B1".to_string(),
958        };
959        let mut values = vec![Some(0.0); 11];
960        values[0] = Some(1.0); // B1
961        values[9] = Some(10.0); // B10
962        values[10] = Some(11.0); // B11 (> B1) => choose B10 = 10
963
964        let result = func.apply(&values).expect("conditional should evaluate");
965        assert_eq!(result, Some(10.0));
966    }
967
968    #[test]
969    fn test_pixel_function_ndvi() {
970        let func = PixelFunction::Ndvi;
971
972        // Standard NDVI calculation
973        let values = vec![Some(0.1), Some(0.5)]; // Red, NIR
974        let result = func.apply(&values);
975        assert!(result.is_ok());
976        // NDVI = (0.5 - 0.1) / (0.5 + 0.1) = 0.4 / 0.6 = 0.666...
977        assert!((result.ok().flatten().expect("Should have value") - 0.666666).abs() < 0.001);
978
979        // With NoData
980        let values_nodata = vec![Some(0.1), None];
981        let result = func.apply(&values_nodata);
982        assert!(result.is_ok());
983        assert_eq!(result.ok().flatten(), None);
984
985        // Zero sum (edge case)
986        let values_zero = vec![Some(0.5), Some(-0.5)];
987        let result = func.apply(&values_zero);
988        assert!(result.is_ok());
989        assert_eq!(result.ok().flatten(), None); // Should return None to avoid division by zero
990    }
991
992    #[test]
993    fn test_pixel_function_evi() {
994        let func = PixelFunction::Evi;
995
996        // Standard EVI calculation
997        let values = vec![Some(0.1), Some(0.5), Some(0.05)]; // Red, NIR, Blue
998        let result = func.apply(&values);
999        assert!(result.is_ok());
1000        // EVI = 2.5 * (0.5 - 0.1) / (0.5 + 6*0.1 - 7.5*0.05 + 1)
1001        // = 2.5 * 0.4 / (0.5 + 0.6 - 0.375 + 1)
1002        // = 1.0 / 1.725 = 0.5797...
1003        let expected = 1.0 / 1.725;
1004        assert!((result.ok().flatten().expect("Should have value") - expected).abs() < 0.001);
1005    }
1006
1007    #[test]
1008    fn test_pixel_function_ndwi() {
1009        let func = PixelFunction::Ndwi;
1010
1011        // Standard NDWI calculation
1012        let values = vec![Some(0.3), Some(0.2)]; // Green, NIR
1013        let result = func.apply(&values);
1014        assert!(result.is_ok());
1015        // NDWI = (0.3 - 0.2) / (0.3 + 0.2) = 0.1 / 0.5 = 0.2
1016        assert!((result.ok().flatten().expect("Should have value") - 0.2).abs() < 0.001);
1017    }
1018
1019    #[test]
1020    fn test_pixel_function_band_math() {
1021        let func = PixelFunction::BandMath {
1022            expression: "(B1 + B2) / 2".to_string(),
1023        };
1024
1025        let values = vec![Some(10.0), Some(20.0)];
1026        let result = func.apply(&values);
1027        assert!(result.is_ok());
1028        assert_eq!(result.ok().flatten(), Some(15.0));
1029
1030        // Test with sqrt
1031        let func_sqrt = PixelFunction::BandMath {
1032            expression: "sqrt(B1)".to_string(),
1033        };
1034        let values_sqrt = vec![Some(16.0)];
1035        let result = func_sqrt.apply(&values_sqrt);
1036        assert!(result.is_ok());
1037        assert_eq!(result.ok().flatten(), Some(4.0));
1038
1039        // Test with abs
1040        let func_abs = PixelFunction::BandMath {
1041            expression: "abs(B1)".to_string(),
1042        };
1043        let values_abs = vec![Some(-5.0)];
1044        let result = func_abs.apply(&values_abs);
1045        assert!(result.is_ok());
1046        assert_eq!(result.ok().flatten(), Some(5.0));
1047    }
1048
1049    #[test]
1050    fn test_pixel_function_lookup_table_nearest() {
1051        let func = PixelFunction::LookupTable {
1052            table: vec![(0.0, 10.0), (0.5, 20.0), (1.0, 30.0)],
1053            interpolation: "nearest".to_string(),
1054        };
1055
1056        // Exact match
1057        let values = vec![Some(0.5)];
1058        let result = func.apply(&values);
1059        assert!(result.is_ok());
1060        assert_eq!(result.ok().flatten(), Some(20.0));
1061
1062        // Nearest neighbor
1063        let values = vec![Some(0.7)];
1064        let result = func.apply(&values);
1065        assert!(result.is_ok());
1066        assert_eq!(result.ok().flatten(), Some(20.0)); // Closest to 0.5
1067    }
1068
1069    #[test]
1070    fn test_pixel_function_lookup_table_linear() {
1071        let func = PixelFunction::LookupTable {
1072            table: vec![(0.0, 10.0), (1.0, 30.0)],
1073            interpolation: "linear".to_string(),
1074        };
1075
1076        // Interpolated value
1077        let values = vec![Some(0.5)];
1078        let result = func.apply(&values);
1079        assert!(result.is_ok());
1080        assert_eq!(result.ok().flatten(), Some(20.0)); // Linear interpolation: 10 + 0.5 * (30-10)
1081
1082        // Edge case: below range
1083        let values = vec![Some(-1.0)];
1084        let result = func.apply(&values);
1085        assert!(result.is_ok());
1086        assert_eq!(result.ok().flatten(), Some(10.0));
1087
1088        // Edge case: above range
1089        let values = vec![Some(2.0)];
1090        let result = func.apply(&values);
1091        assert!(result.is_ok());
1092        assert_eq!(result.ok().flatten(), Some(30.0));
1093    }
1094
1095    #[test]
1096    fn test_pixel_function_conditional() {
1097        let func = PixelFunction::Conditional {
1098            condition: "B1 > 0.5".to_string(),
1099            value_if_true: "B1 * 2".to_string(),
1100            value_if_false: "B1 / 2".to_string(),
1101        };
1102
1103        // True case
1104        let values_true = vec![Some(0.8)];
1105        let result = func.apply(&values_true);
1106        assert!(result.is_ok());
1107        assert_eq!(result.ok().flatten(), Some(1.6));
1108
1109        // False case
1110        let values_false = vec![Some(0.3)];
1111        let result = func.apply(&values_false);
1112        assert!(result.is_ok());
1113        assert_eq!(result.ok().flatten(), Some(0.15));
1114    }
1115
1116    #[test]
1117    fn test_pixel_function_multiply() {
1118        let func = PixelFunction::Multiply;
1119
1120        let values = vec![Some(2.0), Some(3.0), Some(4.0)];
1121        let result = func.apply(&values);
1122        assert!(result.is_ok());
1123        assert_eq!(result.ok().flatten(), Some(24.0));
1124
1125        // With NoData
1126        let values_nodata = vec![Some(2.0), None, Some(4.0)];
1127        let result = func.apply(&values_nodata);
1128        assert!(result.is_ok());
1129        assert_eq!(result.ok().flatten(), Some(8.0)); // Only multiplies valid values
1130    }
1131
1132    #[test]
1133    fn test_pixel_function_divide() {
1134        let func = PixelFunction::Divide;
1135
1136        let values = vec![Some(10.0), Some(2.0)];
1137        let result = func.apply(&values);
1138        assert!(result.is_ok());
1139        assert_eq!(result.ok().flatten(), Some(5.0));
1140
1141        // Division by zero
1142        let values_zero = vec![Some(10.0), Some(0.0)];
1143        let result = func.apply(&values_zero);
1144        assert!(result.is_ok());
1145        assert_eq!(result.ok().flatten(), None);
1146    }
1147
1148    #[test]
1149    fn test_pixel_function_square_root() {
1150        let func = PixelFunction::SquareRoot;
1151
1152        let values = vec![Some(25.0)];
1153        let result = func.apply(&values);
1154        assert!(result.is_ok());
1155        assert_eq!(result.ok().flatten(), Some(5.0));
1156
1157        // Negative value
1158        let values_neg = vec![Some(-4.0)];
1159        let result = func.apply(&values_neg);
1160        assert!(result.is_ok());
1161        assert_eq!(result.ok().flatten(), None);
1162    }
1163
1164    #[test]
1165    fn test_pixel_function_absolute() {
1166        let func = PixelFunction::Absolute;
1167
1168        let values_pos = vec![Some(5.0)];
1169        let result = func.apply(&values_pos);
1170        assert!(result.is_ok());
1171        assert_eq!(result.ok().flatten(), Some(5.0));
1172
1173        let values_neg = vec![Some(-5.0)];
1174        let result = func.apply(&values_neg);
1175        assert!(result.is_ok());
1176        assert_eq!(result.ok().flatten(), Some(5.0));
1177    }
1178
1179    #[test]
1180    fn test_pixel_function_validation() {
1181        // NDVI validation
1182        let ndvi_func = PixelFunction::Ndvi;
1183        let source = VrtSource::new(SourceFilename::absolute("/test.tif"), 1);
1184        let sources_valid = vec![source.clone(), source.clone()];
1185        assert!(ndvi_func.validate(&sources_valid).is_ok());
1186
1187        let sources_invalid = vec![source.clone()];
1188        assert!(ndvi_func.validate(&sources_invalid).is_err());
1189
1190        // BandMath validation
1191        let math_func = PixelFunction::BandMath {
1192            expression: "B1 + B2".to_string(),
1193        };
1194        assert!(math_func.validate(&sources_valid).is_ok());
1195
1196        let empty_expr = PixelFunction::BandMath {
1197            expression: "".to_string(),
1198        };
1199        assert!(empty_expr.validate(&sources_valid).is_err());
1200
1201        // LookupTable validation
1202        let lut_func = PixelFunction::LookupTable {
1203            table: vec![(0.0, 10.0)],
1204            interpolation: "linear".to_string(),
1205        };
1206        assert!(lut_func.validate(&sources_valid).is_ok());
1207
1208        let empty_lut = PixelFunction::LookupTable {
1209            table: vec![],
1210            interpolation: "linear".to_string(),
1211        };
1212        assert!(empty_lut.validate(&sources_valid).is_err());
1213    }
1214}
oxigdal_vrt/band.rs

oxigdal_vrt/
band.rs
