7006568925

Committed 27 Nov 2023 02:07PM UTC coverage: 89.651% (+0.2%) from 89.498%

Build # 7006568925

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #888 from geo-engine/raster_stacks

raster stacking

Run Details

4032 of 4274 new or added lines in 107 files covered. (94.34%)

12 existing lines in 8 files now uncovered.

113020 of 126066 relevant lines covered (89.65%)

59901.79 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

84.29

/operators/src/processing/rgb.rs

use crate::{
    adapters::{QueryWrapper, RasterArrayTimeAdapter},
    engine::{
        BoxRasterQueryProcessor, CanonicOperatorName, ExecutionContext, InitializedRasterOperator,
        InitializedSources, Operator, OperatorData, OperatorName, QueryContext, QueryProcessor,
        RasterBandDescriptors, RasterOperator, RasterQueryProcessor, RasterResultDescriptor,
        TypedRasterQueryProcessor, WorkflowOperatorPath,
    },
    util::Result,
};
use async_trait::async_trait;
use futures::{stream::BoxStream, try_join, StreamExt};
use geoengine_datatypes::{
    dataset::NamedData,
    primitives::{
        partitions_extent, time_interval_extent, BandSelection, RasterQueryRectangle,
        SpatialPartition2D, SpatialResolution,
    },
    raster::{
        FromIndexFn, GridIndexAccess, GridOrEmpty, GridShapeAccess, RasterDataType, RasterTile2D,
    },
    spatial_reference::SpatialReferenceOption,
};
use num_traits::AsPrimitive;
use serde::{Deserialize, Serialize};
use snafu::{ensure, Snafu};

/// The `Rgb` operator combines three raster sources into a single raster output.
/// The output it of type `U32` and contains the red, green and blue values as the first, second and third byte.
/// The forth byte (alpha) is always 255.
pub type Rgb = Operator<RgbParams, RgbSources>;

/// Parameters for the `Rgb` operator.
#[derive(Debug, Serialize, Deserialize, PartialEq, Clone, Copy)]
#[serde(rename_all = "camelCase")]
pub struct RgbParams {
    /// The minimum value for the red channel.
    pub red_min: f64,
    /// The maximum value for the red channel.
    pub red_max: f64,
    /// A scaling factor for the red channel between 0 and 1.
    #[serde(default = "num_traits::One::one")]
    pub red_scale: f64,

    /// The minimum value for the red channel.
    pub green_min: f64,
    /// The maximum value for the red channel.
    pub green_max: f64,
    /// A scaling factor for the green channel between 0 and 1.
    #[serde(default = "num_traits::One::one")]
    pub green_scale: f64,

    /// The minimum value for the red channel.
    pub blue_min: f64,
    /// The maximum value for the red channel.
    pub blue_max: f64,
    /// A scaling factor for the blue channel between 0 and 1.
    #[serde(default = "num_traits::One::one")]
    pub blue_scale: f64,
}

impl RgbParams {
    fn check_valid_user_input(&self) -> Result<()> {
        ensure!(
            self.red_min < self.red_max,
            error::MinGreaterThanMax {
                color: "red",
                min: self.red_min,
                max: self.red_max,
            }
        );
        ensure!(
            self.green_min < self.green_max,
            error::MinGreaterThanMax {
                color: "green",
                min: self.green_min,
                max: self.green_max,
            }
        );
        ensure!(
            self.blue_min < self.blue_max,
            error::MinGreaterThanMax {
                color: "blue",
                min: self.blue_min,
                max: self.blue_max,
            }
        );

        ensure!(
            (0.0..=1.0).contains(&self.red_scale)
                && (0.0..=1.0).contains(&self.green_scale)
                && (0.0..=1.0).contains(&self.blue_scale),
            error::RgbScaleOutOfBounds {
                red: self.red_scale,
                green: self.green_scale,
                blue: self.blue_scale
            }
        );

        Ok(())
    }
}

#[derive(Debug, Clone, Serialize, Deserialize)]
pub struct RgbSources {
    red: Box<dyn RasterOperator>,
    green: Box<dyn RasterOperator>,
    blue: Box<dyn RasterOperator>,
}

impl OperatorData for RgbSources {
    fn data_names_collect(&self, data_names: &mut Vec<NamedData>) {
        for source in self.iter() {
            source.data_names_collect(data_names);
        }
    }
}

impl RgbSources {
    pub fn new(
        red: Box<dyn RasterOperator>,
        green: Box<dyn RasterOperator>,
        blue: Box<dyn RasterOperator>,
    ) -> Self {
        Self { red, green, blue }
    }

    fn iter(&self) -> impl Iterator<Item = &Box<dyn RasterOperator>> {
        [&self.red, &self.green, &self.blue].into_iter()
    }
}

#[async_trait]
impl InitializedSources<RgbInitializedSources> for RgbSources {
    async fn initialize_sources(
        self,
        path: WorkflowOperatorPath,
        context: &dyn ExecutionContext,
    ) -> Result<RgbInitializedSources> {
        let (red, green, blue) = try_join!(
            self.red.initialize(path.clone_and_append(0), context),
            self.green.initialize(path.clone_and_append(1), context),
            self.blue.initialize(path.clone_and_append(2), context),
        )?;

        Ok(RgbInitializedSources { red, green, blue })
    }
}

#[typetag::serde]
#[async_trait]
impl RasterOperator for Rgb {
    async fn _initialize(
        self: Box<Self>,
        path: WorkflowOperatorPath,
        context: &dyn crate::engine::ExecutionContext,
    ) -> Result<Box<dyn InitializedRasterOperator>> {
        self.params.check_valid_user_input()?;

        let name = CanonicOperatorName::from(&self);

        let sources = self.sources.initialize_sources(path, context).await?;

        // TODO: implement multi-band functionality and remove this check
        ensure!(
            sources
                .iter()
                .all(|r| r.result_descriptor().bands.len() == 1),
            crate::error::OperatorDoesNotSupportMultiBandsSourcesYet {
                operator: Rgb::TYPE_NAME
            }
        );

        let spatial_reference = sources.red.result_descriptor().spatial_reference;

        ensure!(
            sources
                .iter()
                .skip(1)
                .all(|rd| rd.result_descriptor().spatial_reference == spatial_reference),
            error::DifferentSpatialReferences {
                red: sources.red.result_descriptor().spatial_reference,
                green: sources.green.result_descriptor().spatial_reference,
                blue: sources.blue.result_descriptor().spatial_reference,
            }
        );

        let time =
            time_interval_extent(sources.iter().map(|source| source.result_descriptor().time));
        let bbox = partitions_extent(sources.iter().map(|source| source.result_descriptor().bbox));

        let resolution = sources
            .iter()
            .map(|source| source.result_descriptor().resolution)
            .reduce(|a, b| match (a, b) {
                (Some(a), Some(b)) => {
                    Some(SpatialResolution::new_unchecked(a.x.min(b.x), a.y.min(b.y)))
                }
                // we can only compute the minimum resolution if all sources have a resolution
                _ => None,
            })
            .flatten();

        let result_descriptor = RasterResultDescriptor {
            data_type: RasterDataType::U32,
            spatial_reference,
            time,
            bbox,
            resolution,
            bands: RasterBandDescriptors::new_single_band(),
        };

        let initialized_operator = InitializedRgb {
            name,
            result_descriptor,
            sources,
            params: self.params,
        };

        Ok(initialized_operator.boxed())
    }

    span_fn!(Rgb);
}

impl OperatorName for Rgb {
    const TYPE_NAME: &'static str = "Rgb";
}

pub struct InitializedRgb {
    name: CanonicOperatorName,
    result_descriptor: RasterResultDescriptor,
    params: RgbParams,
    sources: RgbInitializedSources,
}

pub struct RgbInitializedSources {
    red: Box<dyn InitializedRasterOperator>,
    green: Box<dyn InitializedRasterOperator>,
    blue: Box<dyn InitializedRasterOperator>,
}

impl RgbInitializedSources {
    fn iter(&self) -> impl Iterator<Item = &Box<dyn InitializedRasterOperator>> {
        [&self.red, &self.green, &self.blue].into_iter()
    }
}

impl InitializedRasterOperator for InitializedRgb {
    fn query_processor(&self) -> Result<TypedRasterQueryProcessor> {
        Ok(RgbQueryProcessor::new(
            self.sources.red.query_processor()?.into_f64(),
            self.sources.green.query_processor()?.into_f64(),
            self.sources.blue.query_processor()?.into_f64(),
            self.params,
        )
        .boxed()
        .into())
    }

    fn result_descriptor(&self) -> &RasterResultDescriptor {
        &self.result_descriptor
    }

    fn canonic_name(&self) -> CanonicOperatorName {
        self.name.clone()
    }
}

pub struct RgbQueryProcessor {
    red: BoxRasterQueryProcessor<f64>,
    green: BoxRasterQueryProcessor<f64>,
    blue: BoxRasterQueryProcessor<f64>,
    params: RgbParams,
}

impl RgbQueryProcessor {
    pub fn new(
        red: BoxRasterQueryProcessor<f64>,
        green: BoxRasterQueryProcessor<f64>,
        blue: BoxRasterQueryProcessor<f64>,
        params: RgbParams,
    ) -> Self {
        Self {
            red,
            green,
            blue,
            params,
        }
    }
}

#[async_trait]
impl QueryProcessor for RgbQueryProcessor {
    type Output = RasterTile2D<u32>;
    type SpatialBounds = SpatialPartition2D;
    type Selection = BandSelection;

    async fn _query<'a>(
        &'a self,
        query: RasterQueryRectangle,
        ctx: &'a dyn QueryContext,
    ) -> Result<BoxStream<'a, Result<Self::Output>>> {
        let red = QueryWrapper { p: &self.red, ctx };
        let green = QueryWrapper {
            p: &self.green,
            ctx,
        };
        let blue = QueryWrapper { p: &self.blue, ctx };

        let params = self.params;

        let stream = RasterArrayTimeAdapter::new([red, green, blue], query)
            .map(move |tiles| Ok(compute_tile(tiles?, &params)));

        Ok(Box::pin(stream))
    }
}

fn compute_tile(
    [red, green, blue]: [RasterTile2D<f64>; 3],
    params: &RgbParams,
) -> RasterTile2D<u32> {
    fn fit_to_interval_0_255(value: f64, min: f64, max: f64, scale: f64) -> u32 {
        let mut result = value - min; // shift towards zero
        result /= max - min; // normalize to [0, 1]
        result *= scale; // after scaling with scale ∈ [0, 1], value stays in [0, 1]
        result = (255. * result).round().clamp(0., 255.); // bring value to integer range [0, 255]
        result.as_()
    }

    let map_fn = |lin_idx: usize| -> Option<u32> {
        let (Some(red_value), Some(green_value), Some(blue_value)) = (
            red.get_at_grid_index_unchecked(lin_idx),
            green.get_at_grid_index_unchecked(lin_idx),
            blue.get_at_grid_index_unchecked(lin_idx),
        ) else {
            return None;
        };

        let red =
            fit_to_interval_0_255(red_value, params.red_min, params.red_max, params.red_scale);
        let green = fit_to_interval_0_255(
            green_value,
            params.green_min,
            params.green_max,
            params.green_scale,
        );
        let blue = fit_to_interval_0_255(
            blue_value,
            params.blue_min,
            params.blue_max,
            params.blue_scale,
        );
        let alpha: u32 = 255;

        let rgba = (red << 24) | (green << 16) | (blue << 8) | (alpha);

        Some(rgba)
    };

    // all tiles have the same shape, time, position, etc.
    // so we use red for reference

    let grid_shape = red.grid_shape();
    // TODO: check if parallelism brings any speed-up – might no be faster since we only do simple arithmetic
    let out_grid = GridOrEmpty::from_index_fn(&grid_shape, map_fn);

    RasterTile2D::new(
        red.time,
        red.tile_position,
        0, // TODO
        red.global_geo_transform,
        out_grid,
        red.cache_hint
            .merged(&green.cache_hint)
            .merged(&blue.cache_hint),
    )
}

#[derive(Debug, Snafu)]
#[snafu(visibility(pub(crate)), context(suffix(false)), module(error))]
pub enum RgbOperatorError {
    #[snafu(display(
        "The scale values for r, g and b must be in the range [0, 1]. Got: red: {red}, green: {green}, blue: {blue}",
    ))]
    RgbScaleOutOfBounds { red: f64, green: f64, blue: f64 },

    #[snafu(display(
        "Min value for `{color}` must be smaller than max value. Got: min: {min}, max: {max}",
    ))]
    MinGreaterThanMax {
        color: &'static str,
        min: f64,
        max: f64,
    },

    #[snafu(display(
        "The sources must have the same spatial reference. Got: red: {red}, green: {green}, blue: {blue}",
    ))]
    DifferentSpatialReferences {
        red: SpatialReferenceOption,
        green: SpatialReferenceOption,
        blue: SpatialReferenceOption,
    },
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::engine::{MockExecutionContext, MockQueryContext, QueryProcessor};
    use crate::mock::{MockRasterSource, MockRasterSourceParams};
    use futures::StreamExt;
    use geoengine_datatypes::operations::image::{Colorizer, RgbaColor};
    use geoengine_datatypes::primitives::{
        CacheHint, RasterQueryRectangle, SpatialPartition2D, SpatialResolution, TimeInterval,
    };
    use geoengine_datatypes::raster::{
        Grid2D, GridOrEmpty, MapElements, MaskedGrid2D, RasterTile2D, TileInformation,
        TilingSpecification,
    };
    use geoengine_datatypes::spatial_reference::SpatialReference;
    use geoengine_datatypes::util::test::TestDefault;
    use itertools::Itertools;

    #[test]
    fn deserialize_params() {
        assert_eq!(
            serde_json::from_value::<RgbParams>(serde_json::json!({
                "redMin": 0.,
                "redMax": 1.,
                "redScale": 1.,
                "greenMin": -5.,
                "greenMax": 6.,
                "greenScale": 0.5,
                "blueMin": 0.,
                "blueMax": 10.,
                "blueScale": 0.,
            }))
            .unwrap(),
            RgbParams {
                red_min: 0.,
                red_max: 1.,
                red_scale: 1.,
                green_min: -5.,
                green_max: 6.,
                green_scale: 0.5,
                blue_min: 0.,
                blue_max: 10.,
                blue_scale: 0.,
            }
        );
    }

    #[test]
    fn serialize_params() {
        assert_eq!(
            serde_json::json!({
                "redMin": 0.,
                "redMax": 1.,
                "redScale": 1.,
                "greenMin": -5.,
                "greenMax": 6.,
                "greenScale": 0.5,
                "blueMin": 0.,
                "blueMax": 10.,
                "blueScale": 0.,
            }),
            serde_json::to_value(RgbParams {
                red_min: 0.,
                red_max: 1.,
                red_scale: 1.,
                green_min: -5.,
                green_max: 6.,
                green_scale: 0.5,
                blue_min: 0.,
                blue_max: 10.,
                blue_scale: 0.,
            })
            .unwrap()
        );
    }

    #[tokio::test]
    async fn computation() {
        let tile_size_in_pixels = [3, 2].into();
        let tiling_specification = TilingSpecification {
            origin_coordinate: [0.0, 0.0].into(),
            tile_size_in_pixels,
        };

        let ctx = MockExecutionContext::new_with_tiling_spec(tiling_specification);

        let o = Rgb {
            params: RgbParams {
                red_min: 1.,
                red_max: 6.,
                red_scale: 1.,
                green_min: 1.,
                green_max: 6.,
                green_scale: 0.5,
                blue_min: 1.,
                blue_max: 6.,
                blue_scale: 0.1,
            },
            sources: RgbSources {
                red: make_raster(None),
                green: make_raster(Some(3)),
                blue: make_raster(Some(4)),
            },
        }
        .boxed()
        .initialize(WorkflowOperatorPath::initialize_root(), &ctx)
        .await
        .unwrap();

        let processor = o.query_processor().unwrap().get_u32().unwrap();

        let ctx = MockQueryContext::new(1.into());
        let result_stream = processor
            .query(
                RasterQueryRectangle {
                    spatial_bounds: SpatialPartition2D::new_unchecked(
                        (0., 3.).into(),
                        (2., 0.).into(),
                    ),
                    time_interval: Default::default(),
                    spatial_resolution: SpatialResolution::one(),
                    attributes: BandSelection::first(),
                },
                &ctx,
            )
            .await
            .unwrap();

        let result: Vec<Result<RasterTile2D<u32>>> = result_stream.collect().await;

        assert_eq!(result.len(), 1);

        let first_result = result[0].as_ref().unwrap();

        assert!(!first_result.is_empty());

        let grid = match &first_result.grid_array {
            GridOrEmpty::Grid(g) => g,
            GridOrEmpty::Empty(_) => panic!(),
        };

        let res: Vec<Option<u32>> = grid.masked_element_deref_iterator().collect();

        let expected: [Option<u32>; 6] = [
            Some(0x00_00_00_FF),
            Some(0x33_1A_05_FF),
            None, // 3 is NODATA at green
            None, // 4 is NODATA at blue
            Some(0xCC_66_14_FF),
            Some(0xFF_80_1A_FF),
        ];

        assert_eq!(
            res,
            expected,
            "values must be equal:\n {actual} !=\n {expected}",
            actual = res
                .iter()
                .map(|v| v.map_or("_".to_string(), |v| format!("{v:X}")))
                .join(", "),
            expected = expected
                .iter()
                .map(|v| v.map_or("_".to_string(), |v| format!("{v:X}")))
                .join(", "),
        );

        let colorizer = Colorizer::rgba();
        let color_mapper = colorizer.create_color_mapper();

        let colors: Vec<RgbaColor> = res
            .iter()
            .map(|v| v.map_or(RgbaColor::transparent(), |v| color_mapper.call(v)))
            .collect();

        assert_eq!(
            colors,
            [
                RgbaColor::new(0x00, 0x00, 0x00, 0xFF),
                RgbaColor::new(0x33, 0x1A, 0x05, 0xFF),
                RgbaColor::transparent(), // 3 is NODATA at green
                RgbaColor::transparent(), // 4 is NODATA at blue
                RgbaColor::new(0xCC, 0x66, 0x14, 0xFF),
                RgbaColor::new(0xFF, 0x80, 0x1A, 0xFF),
            ]
        );
    }

    fn make_raster(no_data_value: Option<i8>) -> Box<dyn RasterOperator> {
        let raster = Grid2D::<i8>::new([3, 2].into(), vec![1, 2, 3, 4, 5, 6]).unwrap();

        let real_raster = if let Some(no_data_value) = no_data_value {
            MaskedGrid2D::from(raster)
                .map_elements(|e: Option<i8>| e.filter(|v| *v != no_data_value))
                .into()
        } else {
            GridOrEmpty::from(raster)
        };

        let raster_tile = RasterTile2D::new_with_tile_info(
            TimeInterval::default(),
            TileInformation {
                global_tile_position: [-1, 0].into(),
                tile_size_in_pixels: [3, 2].into(),
                global_geo_transform: TestDefault::test_default(),
            },
            0,
            real_raster,
            CacheHint::no_cache(),
        );

        MockRasterSource {
            params: MockRasterSourceParams {
                data: vec![raster_tile],
                result_descriptor: RasterResultDescriptor {
                    data_type: RasterDataType::I8,
                    spatial_reference: SpatialReference::epsg_4326().into(),
                    time: None,
                    bbox: None,
                    resolution: None,
                    bands: RasterBandDescriptors::new_single_band(),
                },
            },
        }
        .boxed()
    }
}

1	use crate::{
2	adapters::{QueryWrapper, RasterArrayTimeAdapter},
3	engine::{
4	BoxRasterQueryProcessor, CanonicOperatorName, ExecutionContext, InitializedRasterOperator,
5	InitializedSources, Operator, OperatorData, OperatorName, QueryContext, QueryProcessor,
6	RasterBandDescriptors, RasterOperator, RasterQueryProcessor, RasterResultDescriptor,
7	TypedRasterQueryProcessor, WorkflowOperatorPath,
8	},
9	util::Result,
10	};
11	use async_trait::async_trait;
12	use futures::{stream::BoxStream, try_join, StreamExt};
13	use geoengine_datatypes::{
14	dataset::NamedData,
15	primitives::{
16	partitions_extent, time_interval_extent, BandSelection, RasterQueryRectangle,
17	SpatialPartition2D, SpatialResolution,
18	},
19	raster::{
20	FromIndexFn, GridIndexAccess, GridOrEmpty, GridShapeAccess, RasterDataType, RasterTile2D,
21	},
22	spatial_reference::SpatialReferenceOption,
23	};
24	use num_traits::AsPrimitive;
25	use serde::{Deserialize, Serialize};
26	use snafu::{ensure, Snafu};
27
28	/// The `Rgb` operator combines three raster sources into a single raster output.
29	/// The output it of type `U32` and contains the red, green and blue values as the first, second and third byte.
30	/// The forth byte (alpha) is always 255.
31	pub type Rgb = Operator<RgbParams, RgbSources>;
32
33	/// Parameters for the `Rgb` operator.
34	#[derive(Debug, Serialize, Deserialize, PartialEq, Clone, Copy)]	19✔
35	#[serde(rename_all = "camelCase")]
36	pub struct RgbParams {
37	/// The minimum value for the red channel.
38	pub red_min: f64,
39	/// The maximum value for the red channel.
40	pub red_max: f64,
41	/// A scaling factor for the red channel between 0 and 1.
42	#[serde(default = "num_traits::One::one")]
43	pub red_scale: f64,
44
45	/// The minimum value for the red channel.
46	pub green_min: f64,
47	/// The maximum value for the red channel.
48	pub green_max: f64,
49	/// A scaling factor for the green channel between 0 and 1.
50	#[serde(default = "num_traits::One::one")]
51	pub green_scale: f64,
52
53	/// The minimum value for the red channel.
54	pub blue_min: f64,
55	/// The maximum value for the red channel.
56	pub blue_max: f64,
57	/// A scaling factor for the blue channel between 0 and 1.
58	#[serde(default = "num_traits::One::one")]
59	pub blue_scale: f64,
60	}
61
62	impl RgbParams {
63	fn check_valid_user_input(&self) -> Result<()> {	1✔
64	ensure!(	1✔
65	self.red_min < self.red_max,	1✔
66	error::MinGreaterThanMax {	×
67	color: "red",	×
68	min: self.red_min,	×
69	max: self.red_max,	×
70	}	×
71	);
72	ensure!(	1✔
73	self.green_min < self.green_max,	1✔
74	error::MinGreaterThanMax {	×
75	color: "green",	×
76	min: self.green_min,	×
77	max: self.green_max,	×
78	}	×
79	);
80	ensure!(	1✔
81	self.blue_min < self.blue_max,	1✔
82	error::MinGreaterThanMax {	×
83	color: "blue",	×
84	min: self.blue_min,	×
85	max: self.blue_max,	×
86	}	×
87	);
88
89	ensure!(	1✔
90	(0.0..=1.0).contains(&self.red_scale)	1✔
91	&& (0.0..=1.0).contains(&self.green_scale)	1✔
92	&& (0.0..=1.0).contains(&self.blue_scale),	1✔
93	error::RgbScaleOutOfBounds {	×
94	red: self.red_scale,	×
95	green: self.green_scale,	×
96	blue: self.blue_scale	×
97	}	×
98	);
99
100	Ok(())	1✔
101	}	1✔
102	}
103
104	#[derive(Debug, Clone, Serialize, Deserialize)]	1✔
105	pub struct RgbSources {
106	red: Box<dyn RasterOperator>,
107	green: Box<dyn RasterOperator>,
108	blue: Box<dyn RasterOperator>,
109	}
110
111	impl OperatorData for RgbSources {
112	fn data_names_collect(&self, data_names: &mut Vec<NamedData>) {	×
113	for source in self.iter() {	×
114	source.data_names_collect(data_names);	×
115	}	×
116	}	×
117	}
118
119	impl RgbSources {
120	pub fn new(	×
121	red: Box<dyn RasterOperator>,	×
122	green: Box<dyn RasterOperator>,	×
123	blue: Box<dyn RasterOperator>,	×
124	) -> Self {	×
125	Self { red, green, blue }	×
126	}	×
127
128	fn iter(&self) -> impl Iterator<Item = &Box<dyn RasterOperator>> {	×
129	[&self.red, &self.green, &self.blue].into_iter()	×
130	}	×
131	}
132
133	#[async_trait]
134	impl InitializedSources<RgbInitializedSources> for RgbSources {
135	async fn initialize_sources(	1✔
136	self,	1✔
137	path: WorkflowOperatorPath,	1✔
138	context: &dyn ExecutionContext,	1✔
139	) -> Result<RgbInitializedSources> {	1✔
140	let (red, green, blue) = try_join!(	1✔
141	self.red.initialize(path.clone_and_append(0), context),	1✔
142	self.green.initialize(path.clone_and_append(1), context),	1✔
143	self.blue.initialize(path.clone_and_append(2), context),	1✔
144	)?;	1✔
145
146	Ok(RgbInitializedSources { red, green, blue })	1✔
147	}	2✔
148	}
149
150	#[typetag::serde]	×
151	#[async_trait]
152	impl RasterOperator for Rgb {
153	async fn _initialize(	1✔
154	self: Box<Self>,	1✔
155	path: WorkflowOperatorPath,	1✔
156	context: &dyn crate::engine::ExecutionContext,	1✔
157	) -> Result<Box<dyn InitializedRasterOperator>> {	1✔
158	self.params.check_valid_user_input()?;	1✔
159
160	let name = CanonicOperatorName::from(&self);	1✔
161
162	let sources = self.sources.initialize_sources(path, context).await?;	1✔
163
164	// TODO: implement multi-band functionality and remove this check
165	ensure!(	1✔
166	sources	1✔
167	.iter()	1✔
168	.all(\|r\| r.result_descriptor().bands.len() == 1),	3✔
NEW 169	crate::error::OperatorDoesNotSupportMultiBandsSourcesYet {	×
NEW 170	operator: Rgb::TYPE_NAME	×
NEW 171	}	×
172	);
173
174	let spatial_reference = sources.red.result_descriptor().spatial_reference;	1✔
175		1✔
176	ensure!(	1✔
177	sources	1✔
178	.iter()	1✔
179	.skip(1)	1✔
180	.all(\|rd\| rd.result_descriptor().spatial_reference == spatial_reference),	2✔
181	error::DifferentSpatialReferences {	×
182	red: sources.red.result_descriptor().spatial_reference,	×
183	green: sources.green.result_descriptor().spatial_reference,	×
184	blue: sources.blue.result_descriptor().spatial_reference,	×
185	}	×
186	);
187
188	let time =	1✔
189	time_interval_extent(sources.iter().map(\|source\| source.result_descriptor().time));	1✔
190	let bbox = partitions_extent(sources.iter().map(\|source\| source.result_descriptor().bbox));	1✔
191		1✔
192	let resolution = sources	1✔
193	.iter()	1✔
194	.map(\|source\| source.result_descriptor().resolution)	3✔
195	.reduce(\|a, b\| match (a, b) {	2✔
196	(Some(a), Some(b)) => {	×
197	Some(SpatialResolution::new_unchecked(a.x.min(b.x), a.y.min(b.y)))	×
198	}
199	// we can only compute the minimum resolution if all sources have a resolution
200	_ => None,	2✔
201	})	2✔
202	.flatten();	1✔
203		1✔
204	let result_descriptor = RasterResultDescriptor {	1✔
205	data_type: RasterDataType::U32,	1✔
206	spatial_reference,	1✔
207	time,	1✔
208	bbox,	1✔
209	resolution,	1✔
210	bands: RasterBandDescriptors::new_single_band(),	1✔
211	};	1✔
212		1✔
213	let initialized_operator = InitializedRgb {	1✔
214	name,	1✔
215	result_descriptor,	1✔
216	sources,	1✔
217	params: self.params,	1✔
218	};	1✔
219		1✔
220	Ok(initialized_operator.boxed())	1✔
221	}	2✔
222
223	span_fn!(Rgb);	×
224	}
225
226	impl OperatorName for Rgb {
227	const TYPE_NAME: &'static str = "Rgb";
228	}
229
230	pub struct InitializedRgb {
231	name: CanonicOperatorName,
232	result_descriptor: RasterResultDescriptor,
233	params: RgbParams,
234	sources: RgbInitializedSources,
235	}
236
237	pub struct RgbInitializedSources {
238	red: Box<dyn InitializedRasterOperator>,
239	green: Box<dyn InitializedRasterOperator>,
240	blue: Box<dyn InitializedRasterOperator>,
241	}
242
243	impl RgbInitializedSources {
244	fn iter(&self) -> impl Iterator<Item = &Box<dyn InitializedRasterOperator>> {	5✔
245	[&self.red, &self.green, &self.blue].into_iter()	5✔
246	}	5✔
247	}
248
249	impl InitializedRasterOperator for InitializedRgb {
250	fn query_processor(&self) -> Result<TypedRasterQueryProcessor> {	1✔
251	Ok(RgbQueryProcessor::new(	1✔
252	self.sources.red.query_processor()?.into_f64(),	1✔
253	self.sources.green.query_processor()?.into_f64(),	1✔
254	self.sources.blue.query_processor()?.into_f64(),	1✔
255	self.params,	1✔
256	)	1✔
257	.boxed()	1✔
258	.into())	1✔
259	}	1✔
260
261	fn result_descriptor(&self) -> &RasterResultDescriptor {	×
262	&self.result_descriptor	×
263	}	×
264
265	fn canonic_name(&self) -> CanonicOperatorName {	×
266	self.name.clone()	×
267	}	×
268	}
269
270	pub struct RgbQueryProcessor {
271	red: BoxRasterQueryProcessor<f64>,
272	green: BoxRasterQueryProcessor<f64>,
273	blue: BoxRasterQueryProcessor<f64>,
274	params: RgbParams,
275	}
276
277	impl RgbQueryProcessor {
278	pub fn new(	1✔
279	red: BoxRasterQueryProcessor<f64>,	1✔
280	green: BoxRasterQueryProcessor<f64>,	1✔
281	blue: BoxRasterQueryProcessor<f64>,	1✔
282	params: RgbParams,	1✔
283	) -> Self {	1✔
284	Self {	1✔
285	red,	1✔
286	green,	1✔
287	blue,	1✔
288	params,	1✔
289	}	1✔
290	}	1✔
291	}
292
293	#[async_trait]
294	impl QueryProcessor for RgbQueryProcessor {
295	type Output = RasterTile2D<u32>;
296	type SpatialBounds = SpatialPartition2D;
297	type Selection = BandSelection;
298
299	async fn _query<'a>(	1✔
300	&'a self,	1✔
301	query: RasterQueryRectangle,	1✔
302	ctx: &'a dyn QueryContext,	1✔
303	) -> Result<BoxStream<'a, Result<Self::Output>>> {	1✔
304	let red = QueryWrapper { p: &self.red, ctx };	1✔
305	let green = QueryWrapper {	1✔
306	p: &self.green,	1✔
307	ctx,	1✔
308	};	1✔
309	let blue = QueryWrapper { p: &self.blue, ctx };	1✔
310		1✔
311	let params = self.params;	1✔
312		1✔
313	let stream = RasterArrayTimeAdapter::new([red, green, blue], query)	1✔
314	.map(move \|tiles\| Ok(compute_tile(tiles?, &params)));	1✔
315		1✔
316	Ok(Box::pin(stream))	1✔
317	}	1✔
318	}
319
320	fn compute_tile(	1✔
321	[red, green, blue]: [RasterTile2D<f64>; 3],	1✔
322	params: &RgbParams,	1✔
323	) -> RasterTile2D<u32> {	1✔
324	fn fit_to_interval_0_255(value: f64, min: f64, max: f64, scale: f64) -> u32 {	12✔
325	let mut result = value - min; // shift towards zero	12✔
326	result /= max - min; // normalize to [0, 1]	12✔
327	result *= scale; // after scaling with scale ∈ [0, 1], value stays in [0, 1]	12✔
328	result = (255. * result).round().clamp(0., 255.); // bring value to integer range [0, 255]	12✔
329	result.as_()	12✔
330	}	12✔
331		1✔
332	let map_fn = \|lin_idx: usize\| -> Option<u32> {	6✔
333	let (Some(red_value), Some(green_value), Some(blue_value)) = (	4✔
334	red.get_at_grid_index_unchecked(lin_idx),	6✔
335	green.get_at_grid_index_unchecked(lin_idx),	6✔
336	blue.get_at_grid_index_unchecked(lin_idx),	6✔
337	) else {
338	return None;	2✔
339	};
340
341	let red =	4✔
342	fit_to_interval_0_255(red_value, params.red_min, params.red_max, params.red_scale);	4✔
343	let green = fit_to_interval_0_255(	4✔
344	green_value,	4✔
345	params.green_min,	4✔
346	params.green_max,	4✔
347	params.green_scale,	4✔
348	);	4✔
349	let blue = fit_to_interval_0_255(	4✔
350	blue_value,	4✔
351	params.blue_min,	4✔
352	params.blue_max,	4✔
353	params.blue_scale,	4✔
354	);	4✔
355	let alpha: u32 = 255;	4✔
356		4✔
357	let rgba = (red << 24) \| (green << 16) \| (blue << 8) \| (alpha);	4✔
358		4✔
359	Some(rgba)	4✔
360	};	6✔
361		1✔
362	// all tiles have the same shape, time, position, etc.	1✔
363	// so we use red for reference	1✔
364		1✔
365	let grid_shape = red.grid_shape();	1✔
366	// TODO: check if parallelism brings any speed-up – might no be faster since we only do simple arithmetic	1✔
367	let out_grid = GridOrEmpty::from_index_fn(&grid_shape, map_fn);	1✔
368		1✔
369	RasterTile2D::new(	1✔
370	red.time,	1✔
371	red.tile_position,	1✔
372	0, // TODO	1✔
373	red.global_geo_transform,	1✔
374	out_grid,	1✔
375	red.cache_hint	1✔
376	.merged(&green.cache_hint)	1✔
377	.merged(&blue.cache_hint),	1✔
378	)	1✔
379	}	1✔
380
381	#[derive(Debug, Snafu)]	×
382	#[snafu(visibility(pub(crate)), context(suffix(false)), module(error))]
383	pub enum RgbOperatorError {
384	#[snafu(display(
385	"The scale values for r, g and b must be in the range [0, 1]. Got: red: {red}, green: {green}, blue: {blue}",
386	))]
387	RgbScaleOutOfBounds { red: f64, green: f64, blue: f64 },
388
389	#[snafu(display(
390	"Min value for `{color}` must be smaller than max value. Got: min: {min}, max: {max}",
391	))]
392	MinGreaterThanMax {
393	color: &'static str,
394	min: f64,
395	max: f64,
396	},
397
398	#[snafu(display(
399	"The sources must have the same spatial reference. Got: red: {red}, green: {green}, blue: {blue}",
400	))]
401	DifferentSpatialReferences {
402	red: SpatialReferenceOption,
403	green: SpatialReferenceOption,
404	blue: SpatialReferenceOption,
405	},
406	}
407
408	#[cfg(test)]
409	mod tests {
410	use super::*;
411	use crate::engine::{MockExecutionContext, MockQueryContext, QueryProcessor};
412	use crate::mock::{MockRasterSource, MockRasterSourceParams};
413	use futures::StreamExt;
414	use geoengine_datatypes::operations::image::{Colorizer, RgbaColor};
415	use geoengine_datatypes::primitives::{
416	CacheHint, RasterQueryRectangle, SpatialPartition2D, SpatialResolution, TimeInterval,
417	};
418	use geoengine_datatypes::raster::{
419	Grid2D, GridOrEmpty, MapElements, MaskedGrid2D, RasterTile2D, TileInformation,
420	TilingSpecification,
421	};
422	use geoengine_datatypes::spatial_reference::SpatialReference;
423	use geoengine_datatypes::util::test::TestDefault;
424	use itertools::Itertools;
425
426	#[test]	1✔
427	fn deserialize_params() {	1✔
428	assert_eq!(	1✔
429	serde_json::from_value::<RgbParams>(serde_json::json!({	1✔
430	"redMin": 0.,	1✔
431	"redMax": 1.,	1✔
432	"redScale": 1.,	1✔
433	"greenMin": -5.,	1✔
434	"greenMax": 6.,	1✔
435	"greenScale": 0.5,	1✔
436	"blueMin": 0.,	1✔
437	"blueMax": 10.,	1✔
438	"blueScale": 0.,	1✔
439	}))	1✔
440	.unwrap(),	1✔
441	RgbParams {	1✔
442	red_min: 0.,	1✔
443	red_max: 1.,	1✔
444	red_scale: 1.,	1✔
445	green_min: -5.,	1✔
446	green_max: 6.,	1✔
447	green_scale: 0.5,	1✔
448	blue_min: 0.,	1✔
449	blue_max: 10.,	1✔
450	blue_scale: 0.,	1✔
451	}	1✔
452	);	1✔
453	}	1✔
454
455	#[test]	1✔
456	fn serialize_params() {	1✔
457	assert_eq!(	1✔
458	serde_json::json!({	1✔
459	"redMin": 0.,	1✔
460	"redMax": 1.,	1✔
461	"redScale": 1.,	1✔
462	"greenMin": -5.,	1✔
463	"greenMax": 6.,	1✔
464	"greenScale": 0.5,	1✔
465	"blueMin": 0.,	1✔
466	"blueMax": 10.,	1✔
467	"blueScale": 0.,	1✔
468	}),	1✔
469	serde_json::to_value(RgbParams {	1✔
470	red_min: 0.,	1✔
471	red_max: 1.,	1✔
472	red_scale: 1.,	1✔
473	green_min: -5.,	1✔
474	green_max: 6.,	1✔
475	green_scale: 0.5,	1✔
476	blue_min: 0.,	1✔
477	blue_max: 10.,	1✔
478	blue_scale: 0.,	1✔
479	})	1✔
480	.unwrap()	1✔
481	);	1✔
482	}	1✔
483
484	#[tokio::test]	1✔
485	async fn computation() {	1✔
486	let tile_size_in_pixels = [3, 2].into();	1✔
487	let tiling_specification = TilingSpecification {	1✔
488	origin_coordinate: [0.0, 0.0].into(),	1✔
489	tile_size_in_pixels,	1✔
490	};	1✔
491		1✔
492	let ctx = MockExecutionContext::new_with_tiling_spec(tiling_specification);	1✔
493
494	let o = Rgb {	1✔
495	params: RgbParams {	1✔
496	red_min: 1.,	1✔
497	red_max: 6.,	1✔
498	red_scale: 1.,	1✔
499	green_min: 1.,	1✔
500	green_max: 6.,	1✔
501	green_scale: 0.5,	1✔
502	blue_min: 1.,	1✔
503	blue_max: 6.,	1✔
504	blue_scale: 0.1,	1✔
505	},	1✔
506	sources: RgbSources {	1✔
507	red: make_raster(None),	1✔
508	green: make_raster(Some(3)),	1✔
509	blue: make_raster(Some(4)),	1✔
510	},	1✔
511	}	1✔
512	.boxed()	1✔
513	.initialize(WorkflowOperatorPath::initialize_root(), &ctx)	1✔
514	.await	×
515	.unwrap();	1✔
516		1✔
517	let processor = o.query_processor().unwrap().get_u32().unwrap();	1✔
518		1✔
519	let ctx = MockQueryContext::new(1.into());	1✔
520	let result_stream = processor	1✔
521	.query(	1✔
522	RasterQueryRectangle {	1✔
523	spatial_bounds: SpatialPartition2D::new_unchecked(	1✔
524	(0., 3.).into(),	1✔
525	(2., 0.).into(),	1✔
526	),	1✔
527	time_interval: Default::default(),	1✔
528	spatial_resolution: SpatialResolution::one(),	1✔
529	attributes: BandSelection::first(),	1✔
530	},	1✔
531	&ctx,	1✔
532	)	1✔
533	.await	×
534	.unwrap();	1✔
535
536	let result: Vec<Result<RasterTile2D<u32>>> = result_stream.collect().await;	1✔
537
538	assert_eq!(result.len(), 1);	1✔
539
540	let first_result = result[0].as_ref().unwrap();	1✔
541		1✔
542	assert!(!first_result.is_empty());	1✔
543
544	let grid = match &first_result.grid_array {	1✔
545	GridOrEmpty::Grid(g) => g,	1✔
546	GridOrEmpty::Empty(_) => panic!(),	×
547	};
548
549	let res: Vec<Option<u32>> = grid.masked_element_deref_iterator().collect();	1✔
550		1✔
551	let expected: [Option<u32>; 6] = [	1✔
552	Some(0x00_00_00_FF),	1✔
553	Some(0x33_1A_05_FF),	1✔
554	None, // 3 is NODATA at green	1✔
555	None, // 4 is NODATA at blue	1✔
556	Some(0xCC_66_14_FF),	1✔
557	Some(0xFF_80_1A_FF),	1✔
558	];	1✔
559		1✔
560	assert_eq!(	1✔
561	res,
562	expected,
563	"values must be equal:\n {actual} !=\n {expected}",	×
564	actual = res	×
565	.iter()	×
566	.map(\|v\| v.map_or("_".to_string(), \|v\| format!("{v:X}")))	×
567	.join(", "),	×
568	expected = expected	×
569	.iter()	×
570	.map(\|v\| v.map_or("_".to_string(), \|v\| format!("{v:X}")))	×
571	.join(", "),	×
572	);
573
574	let colorizer = Colorizer::rgba();	1✔
575	let color_mapper = colorizer.create_color_mapper();	1✔
576		1✔
577	let colors: Vec<RgbaColor> = res	1✔
578	.iter()	1✔
579	.map(\|v\| v.map_or(RgbaColor::transparent(), \|v\| color_mapper.call(v)))	6✔
580	.collect();	1✔
581		1✔
582	assert_eq!(	1✔
583	colors,	1✔
584	[	1✔
585	RgbaColor::new(0x00, 0x00, 0x00, 0xFF),	1✔
586	RgbaColor::new(0x33, 0x1A, 0x05, 0xFF),	1✔
587	RgbaColor::transparent(), // 3 is NODATA at green	1✔
588	RgbaColor::transparent(), // 4 is NODATA at blue	1✔
589	RgbaColor::new(0xCC, 0x66, 0x14, 0xFF),	1✔
590	RgbaColor::new(0xFF, 0x80, 0x1A, 0xFF),	1✔
591	]	1✔
592	);	1✔
593	}
594
595	fn make_raster(no_data_value: Option<i8>) -> Box<dyn RasterOperator> {	3✔
596	let raster = Grid2D::<i8>::new([3, 2].into(), vec![1, 2, 3, 4, 5, 6]).unwrap();	3✔
597
598	let real_raster = if let Some(no_data_value) = no_data_value {	3✔
599	MaskedGrid2D::from(raster)	2✔
600	.map_elements(\|e: Option<i8>\| e.filter(\|v\| *v != no_data_value))	12✔
601	.into()	2✔
602	} else {
603	GridOrEmpty::from(raster)	1✔
604	};
605
606	let raster_tile = RasterTile2D::new_with_tile_info(	3✔
607	TimeInterval::default(),	3✔
608	TileInformation {	3✔
609	global_tile_position: [-1, 0].into(),	3✔
610	tile_size_in_pixels: [3, 2].into(),	3✔
611	global_geo_transform: TestDefault::test_default(),	3✔
612	},	3✔
613	0,	3✔
614	real_raster,	3✔
615	CacheHint::no_cache(),	3✔
616	);	3✔
617		3✔
618	MockRasterSource {	3✔
619	params: MockRasterSourceParams {	3✔
620	data: vec![raster_tile],	3✔
621	result_descriptor: RasterResultDescriptor {	3✔
622	data_type: RasterDataType::I8,	3✔
623	spatial_reference: SpatialReference::epsg_4326().into(),	3✔
624	time: None,	3✔
625	bbox: None,	3✔
626	resolution: None,	3✔
627	bands: RasterBandDescriptors::new_single_band(),	3✔
628	},	3✔
629	},	3✔
630	}	3✔
631	.boxed()	3✔
632	}	3✔
633	}

geo-engine / geoengine / 7006568925

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous