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

93.83

/operators/src/plot/temporal_raster_mean_plot.rs

use crate::engine::{
    CanonicOperatorName, ExecutionContext, InitializedPlotOperator, InitializedRasterOperator,
    InitializedSources, Operator, OperatorName, PlotOperator, PlotQueryProcessor,
    PlotResultDescriptor, QueryContext, QueryProcessor, RasterQueryProcessor, SingleRasterSource,
    TypedPlotQueryProcessor, WorkflowOperatorPath,
};
use crate::util::math::average_floor;
use crate::util::Result;
use async_trait::async_trait;
use futures::stream::BoxStream;
use futures::StreamExt;
use geoengine_datatypes::plots::{AreaLineChart, Plot, PlotData};
use geoengine_datatypes::primitives::{
    BandSelection, Measurement, PlotQueryRectangle, RasterQueryRectangle, TimeInstance,
    TimeInterval,
};
use geoengine_datatypes::raster::{Pixel, RasterTile2D};
use serde::{Deserialize, Serialize};
use snafu::ensure;
use std::collections::BTreeMap;

pub const MEAN_RASTER_PIXEL_VALUES_OVER_TIME_NAME: &str = "Mean Raster Pixel Values over Time";

/// A plot that shows the mean values of rasters over time as an area plot.
pub type MeanRasterPixelValuesOverTime =
    Operator<MeanRasterPixelValuesOverTimeParams, SingleRasterSource>;

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

/// The parameter spec for `MeanRasterPixelValuesOverTime`
#[derive(Debug, Clone, PartialEq, Eq, Serialize, Deserialize)]
#[serde(rename_all = "camelCase")]
pub struct MeanRasterPixelValuesOverTimeParams {
    /// Where should the x-axis (time) tick be positioned?
    /// At either time start, time end or in the center.
    pub time_position: MeanRasterPixelValuesOverTimePosition,

    /// Whether to fill the area under the curve.
    #[serde(default = "default_true")]
    pub area: bool,
}

const fn default_true() -> bool {
    true
}

#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
#[serde(rename_all = "camelCase")]
pub enum MeanRasterPixelValuesOverTimePosition {
    Start,
    Center,
    End,
}

#[typetag::serde]
#[async_trait]
impl PlotOperator for MeanRasterPixelValuesOverTime {
    async fn _initialize(
        self: Box<Self>,
        path: WorkflowOperatorPath,
        context: &dyn ExecutionContext,
    ) -> Result<Box<dyn InitializedPlotOperator>> {
        let name = CanonicOperatorName::from(&self);

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

        let in_desc = raster.result_descriptor().clone();

        // TODO: implement multi-band functionality and remove this check
        ensure!(
            in_desc.bands.len() == 1,
            crate::error::OperatorDoesNotSupportMultiBandsSourcesYet {
                operator: MeanRasterPixelValuesOverTime::TYPE_NAME
            }
        );

        let initialized_operator = InitializedMeanRasterPixelValuesOverTime {
            name,
            result_descriptor: in_desc.into(),
            raster,
            state: self.params,
        };

        Ok(initialized_operator.boxed())
    }

    span_fn!(MeanRasterPixelValuesOverTime);
}

/// The initialization of `MeanRasterPixelValuesOverTime`
pub struct InitializedMeanRasterPixelValuesOverTime {
    name: CanonicOperatorName,
    result_descriptor: PlotResultDescriptor,
    raster: Box<dyn InitializedRasterOperator>,
    state: MeanRasterPixelValuesOverTimeParams,
}

impl InitializedPlotOperator for InitializedMeanRasterPixelValuesOverTime {
    fn query_processor(&self) -> Result<TypedPlotQueryProcessor> {
        let input_processor = self.raster.query_processor()?;
        let time_position = self.state.time_position;
        let measurement = self.raster.result_descriptor().bands[0].measurement.clone(); // TODO: adjust for multibands
        let draw_area = self.state.area;

        let processor = call_on_generic_raster_processor!(input_processor, raster => {
            MeanRasterPixelValuesOverTimeQueryProcessor { raster, time_position, measurement, draw_area }.boxed()
        });

        Ok(TypedPlotQueryProcessor::JsonVega(processor))
    }

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

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

/// A query processor that calculates the `TemporalRasterMeanPlot` about its input.
pub struct MeanRasterPixelValuesOverTimeQueryProcessor<P: Pixel> {
    raster: Box<dyn RasterQueryProcessor<RasterType = P>>,
    time_position: MeanRasterPixelValuesOverTimePosition,
    measurement: Measurement,
    draw_area: bool,
}

#[async_trait]
impl<P: Pixel> PlotQueryProcessor for MeanRasterPixelValuesOverTimeQueryProcessor<P> {
    type OutputFormat = PlotData;

    fn plot_type(&self) -> &'static str {
        MEAN_RASTER_PIXEL_VALUES_OVER_TIME_NAME
    }

    async fn plot_query<'a>(
        &'a self,
        query: PlotQueryRectangle,
        ctx: &'a dyn QueryContext,
    ) -> Result<Self::OutputFormat> {
        let means = Self::calculate_means(
            self.raster
                .query(
                    RasterQueryRectangle::from_qrect_and_bands(&query, BandSelection::first()),
                    ctx,
                )
                .await?,
            self.time_position,
        )
        .await?;

        let plot = Self::generate_plot(means, self.measurement.clone(), self.draw_area)?;

        let plot_data = plot.to_vega_embeddable(false)?;

        Ok(plot_data)
    }
}

impl<P: Pixel> MeanRasterPixelValuesOverTimeQueryProcessor<P> {
    async fn calculate_means(
        mut tile_stream: BoxStream<'_, Result<RasterTile2D<P>>>,
        position: MeanRasterPixelValuesOverTimePosition,
    ) -> Result<BTreeMap<TimeInstance, MeanCalculator>> {
        let mut means: BTreeMap<TimeInstance, MeanCalculator> = BTreeMap::new();

        while let Some(tile) = tile_stream.next().await {
            let tile = tile?;

            match tile.grid_array {
                geoengine_datatypes::raster::GridOrEmpty::Grid(g) => {
                    let time = Self::time_interval_projection(tile.time, position);
                    let mean = means.entry(time).or_default();
                    mean.add(g.masked_element_deref_iterator());
                }
                geoengine_datatypes::raster::GridOrEmpty::Empty(_) => (),
            }
        }

        Ok(means)
    }

    #[inline]
    fn time_interval_projection(
        time_interval: TimeInterval,
        time_position: MeanRasterPixelValuesOverTimePosition,
    ) -> TimeInstance {
        match time_position {
            MeanRasterPixelValuesOverTimePosition::Start => time_interval.start(),
            MeanRasterPixelValuesOverTimePosition::Center => TimeInstance::from_millis_unchecked(
                average_floor(time_interval.start().inner(), time_interval.end().inner()),
            ),
            MeanRasterPixelValuesOverTimePosition::End => time_interval.end(),
        }
    }

    fn generate_plot(
        means: BTreeMap<TimeInstance, MeanCalculator>,
        measurement: Measurement,
        draw_area: bool,
    ) -> Result<AreaLineChart> {
        let mut timestamps = Vec::with_capacity(means.len());
        let mut values = Vec::with_capacity(means.len());

        for (timestamp, mean_calculator) in means {
            timestamps.push(timestamp);
            values.push(mean_calculator.mean());
        }

        AreaLineChart::new(timestamps, values, measurement, draw_area).map_err(Into::into)
    }
}

struct MeanCalculator {
    mean: f64,
    n: usize,
}

impl Default for MeanCalculator {
    fn default() -> Self {
        Self { mean: 0., n: 0 }
    }
}

impl MeanCalculator {
    #[inline]
    fn add<P: Pixel, I: Iterator<Item = Option<P>>>(&mut self, values: I) {
        values.flatten().for_each(|value| {
            self.add_single_value(value);
        });
    }

    #[inline]
    fn add_single_value<P: Pixel>(&mut self, value: P) {
        let value: f64 = value.as_();

        if value.is_nan() {
            return;
        }

        self.n += 1;
        let delta = value - self.mean;
        self.mean += delta / (self.n as f64);
    }

    #[inline]
    fn mean(&self) -> f64 {
        self.mean
    }
}

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

    use crate::{
        engine::{
            ChunkByteSize, MockExecutionContext, MockQueryContext, RasterBandDescriptors,
            RasterOperator, RasterResultDescriptor,
        },
        source::GdalSource,
    };
    use crate::{
        mock::{MockRasterSource, MockRasterSourceParams},
        source::GdalSourceParameters,
    };
    use geoengine_datatypes::primitives::{
        BoundingBox2D, CacheHint, Measurement, PlotSeriesSelection, SpatialResolution, TimeInterval,
    };
    use geoengine_datatypes::{dataset::NamedData, plots::PlotMetaData, primitives::DateTime};
    use geoengine_datatypes::{raster::TilingSpecification, spatial_reference::SpatialReference};
    use geoengine_datatypes::{
        raster::{Grid2D, RasterDataType, TileInformation},
        util::test::TestDefault,
    };
    use serde_json::{json, Value};

    #[test]
    fn serialization() {
        let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {
            params: MeanRasterPixelValuesOverTimeParams {
                time_position: MeanRasterPixelValuesOverTimePosition::Start,
                area: true,
            },
            sources: SingleRasterSource {
                raster: GdalSource {
                    params: GdalSourceParameters {
                        data: NamedData::with_system_name("test"),
                    },
                }
                .boxed(),
            },
        };

        let serialized = json!({
            "type": "MeanRasterPixelValuesOverTime",
            "params": {
                "timePosition": "start",
                "area": true,
            },
            "sources": {
                "raster": {
                    "type": "GdalSource",
                    "params": {
                        "data": "test"
                    }
                }
            },
            "vectorSources": [],
        })
        .to_string();

        serde_json::from_str::<Box<dyn PlotOperator>>(&serialized).unwrap();

        let deserialized: MeanRasterPixelValuesOverTime =
            serde_json::from_str(&serialized).unwrap();

        assert_eq!(deserialized.params, temporal_raster_mean_plot.params);
    }

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

        let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {
            params: MeanRasterPixelValuesOverTimeParams {
                time_position: MeanRasterPixelValuesOverTimePosition::Center,
                area: true,
            },
            sources: SingleRasterSource {
                raster: generate_mock_raster_source(
                    vec![TimeInterval::new(
                        TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),
                        TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),
                    )
                    .unwrap()],
                    vec![vec![1, 2, 3, 4, 5, 6]],
                ),
            },
        };

        let temporal_raster_mean_plot = temporal_raster_mean_plot
            .boxed()
            .initialize(WorkflowOperatorPath::initialize_root(), &execution_context)
            .await
            .unwrap();

        let processor = temporal_raster_mean_plot
            .query_processor()
            .unwrap()
            .json_vega()
            .unwrap();

        let result = processor
            .plot_query(
                PlotQueryRectangle {
                    spatial_bounds: BoundingBox2D::new((-180., -90.).into(), (180., 90.).into())
                        .unwrap(),
                    time_interval: TimeInterval::default(),
                    spatial_resolution: SpatialResolution::one(),
                    attributes: PlotSeriesSelection::all(),
                },
                &MockQueryContext::new(ChunkByteSize::MIN),
            )
            .await
            .unwrap();

        assert!(matches!(result.metadata, PlotMetaData::None));

        let vega_json: Value = serde_json::from_str(&result.vega_string).unwrap();

        assert_eq!(
            vega_json,
            json!({
                "$schema": "https://vega.github.io/schema/vega-lite/v4.17.0.json",
                "data": {
                    "values": [{
                        "x": "1995-01-01T00:00:00+00:00",
                        "y": 3.5
                    }]
                },
                "description": "Area Plot",
                "encoding": {
                    "x": {
                        "field": "x",
                        "title": "Time",
                        "type": "temporal"
                    },
                    "y": {
                        "field": "y",
                        "title": "",
                        "type": "quantitative"
                    }
                },
                "mark": {
                    "type": "area",
                    "line": true,
                    "point": true
                }
            })
        );
    }

    fn generate_mock_raster_source(
        time_intervals: Vec<TimeInterval>,
        values_vec: Vec<Vec<u8>>,
    ) -> Box<dyn RasterOperator> {
        assert_eq!(time_intervals.len(), values_vec.len());
        assert!(values_vec.iter().all(|v| v.len() == 6));

        let mut tiles = Vec::with_capacity(time_intervals.len());
        for (time_interval, values) in time_intervals.into_iter().zip(values_vec) {
            tiles.push(RasterTile2D::new_with_tile_info(
                time_interval,
                TileInformation {
                    global_geo_transform: TestDefault::test_default(),
                    global_tile_position: [0, 0].into(),
                    tile_size_in_pixels: [3, 2].into(),
                },
                0,
                Grid2D::new([3, 2].into(), values).unwrap().into(),
                CacheHint::default(),
            ));
        }

        MockRasterSource {
            params: MockRasterSourceParams {
                data: tiles,
                result_descriptor: RasterResultDescriptor {
                    data_type: RasterDataType::U8,
                    spatial_reference: SpatialReference::epsg_4326().into(),
                    time: None,
                    bbox: None,
                    resolution: None,
                    bands: RasterBandDescriptors::new_single_band(),
                },
            },
        }
        .boxed()
    }

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

        let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {
            params: MeanRasterPixelValuesOverTimeParams {
                time_position: MeanRasterPixelValuesOverTimePosition::Start,
                area: true,
            },

            sources: SingleRasterSource {
                raster: generate_mock_raster_source(
                    vec![
                        TimeInterval::new(
                            TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),
                            TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),
                        )
                        .unwrap(),
                        TimeInterval::new(
                            TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),
                            TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),
                        )
                        .unwrap(),
                        TimeInterval::new(
                            TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),
                            TimeInstance::from(DateTime::new_utc(2005, 1, 1, 0, 0, 0)),
                        )
                        .unwrap(),
                    ],
                    vec![
                        vec![1, 2, 3, 4, 5, 6],
                        vec![9, 9, 8, 8, 8, 9],
                        vec![3, 4, 5, 6, 7, 8],
                    ],
                ),
            },
        };

        let temporal_raster_mean_plot = temporal_raster_mean_plot
            .boxed()
            .initialize(WorkflowOperatorPath::initialize_root(), &execution_context)
            .await
            .unwrap();

        let processor = temporal_raster_mean_plot
            .query_processor()
            .unwrap()
            .json_vega()
            .unwrap();

        let result = processor
            .plot_query(
                PlotQueryRectangle {
                    spatial_bounds: BoundingBox2D::new((-180., -90.).into(), (180., 90.).into())
                        .unwrap(),
                    time_interval: TimeInterval::default(),
                    spatial_resolution: SpatialResolution::one(),
                    attributes: PlotSeriesSelection::all(),
                },
                &MockQueryContext::new(ChunkByteSize::MIN),
            )
            .await
            .unwrap();

        assert_eq!(
            result,
            AreaLineChart::new(
                vec![
                    TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),
                    TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),
                    TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0))
                ],
                vec![3.5, 8.5, 5.5],
                Measurement::Unitless,
                true,
            )
            .unwrap()
            .to_vega_embeddable(false)
            .unwrap()
        );
    }
}

1	use crate::engine::{
2	CanonicOperatorName, ExecutionContext, InitializedPlotOperator, InitializedRasterOperator,
3	InitializedSources, Operator, OperatorName, PlotOperator, PlotQueryProcessor,
4	PlotResultDescriptor, QueryContext, QueryProcessor, RasterQueryProcessor, SingleRasterSource,
5	TypedPlotQueryProcessor, WorkflowOperatorPath,
6	};
7	use crate::util::math::average_floor;
8	use crate::util::Result;
9	use async_trait::async_trait;
10	use futures::stream::BoxStream;
11	use futures::StreamExt;
12	use geoengine_datatypes::plots::{AreaLineChart, Plot, PlotData};
13	use geoengine_datatypes::primitives::{
14	BandSelection, Measurement, PlotQueryRectangle, RasterQueryRectangle, TimeInstance,
15	TimeInterval,
16	};
17	use geoengine_datatypes::raster::{Pixel, RasterTile2D};
18	use serde::{Deserialize, Serialize};
19	use snafu::ensure;
20	use std::collections::BTreeMap;
21
22	pub const MEAN_RASTER_PIXEL_VALUES_OVER_TIME_NAME: &str = "Mean Raster Pixel Values over Time";
23
24	/// A plot that shows the mean values of rasters over time as an area plot.
25	pub type MeanRasterPixelValuesOverTime =
26	Operator<MeanRasterPixelValuesOverTimeParams, SingleRasterSource>;
27
28	impl OperatorName for MeanRasterPixelValuesOverTime {
29	const TYPE_NAME: &'static str = "MeanRasterPixelValuesOverTime";
30	}
31
32	/// The parameter spec for `MeanRasterPixelValuesOverTime`
33	#[derive(Debug, Clone, PartialEq, Eq, Serialize, Deserialize)]	10✔
34	#[serde(rename_all = "camelCase")]
35	pub struct MeanRasterPixelValuesOverTimeParams {
36	/// Where should the x-axis (time) tick be positioned?
37	/// At either time start, time end or in the center.
38	pub time_position: MeanRasterPixelValuesOverTimePosition,
39
40	/// Whether to fill the area under the curve.
41	#[serde(default = "default_true")]
42	pub area: bool,
43	}
44
45	const fn default_true() -> bool {	×
46	true	×
47	}	×
48
49	#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]	4✔
50	#[serde(rename_all = "camelCase")]
51	pub enum MeanRasterPixelValuesOverTimePosition {
52	Start,
53	Center,
54	End,
55	}
56
57	#[typetag::serde]	1✔
58	#[async_trait]
59	impl PlotOperator for MeanRasterPixelValuesOverTime {
60	async fn _initialize(	3✔
61	self: Box<Self>,	3✔
62	path: WorkflowOperatorPath,	3✔
63	context: &dyn ExecutionContext,	3✔
64	) -> Result<Box<dyn InitializedPlotOperator>> {	3✔
65	let name = CanonicOperatorName::from(&self);	3✔
66
67	let initalized_sources = self.sources.initialize_sources(path, context).await?;	10✔
68	let raster = initalized_sources.raster;	3✔
69		3✔
70	let in_desc = raster.result_descriptor().clone();	3✔
71		3✔
72	// TODO: implement multi-band functionality and remove this check	3✔
73	ensure!(	3✔
74	in_desc.bands.len() == 1,	3✔
NEW 75	crate::error::OperatorDoesNotSupportMultiBandsSourcesYet {	×
NEW 76	operator: MeanRasterPixelValuesOverTime::TYPE_NAME	×
NEW 77	}	×
78	);
79
80	let initialized_operator = InitializedMeanRasterPixelValuesOverTime {	3✔
81	name,	3✔
82	result_descriptor: in_desc.into(),	3✔
83	raster,	3✔
84	state: self.params,	3✔
85	};	3✔
86		3✔
87	Ok(initialized_operator.boxed())	3✔
88	}	6✔
89
90	span_fn!(MeanRasterPixelValuesOverTime);	×
91	}
92
93	/// The initialization of `MeanRasterPixelValuesOverTime`
94	pub struct InitializedMeanRasterPixelValuesOverTime {
95	name: CanonicOperatorName,
96	result_descriptor: PlotResultDescriptor,
97	raster: Box<dyn InitializedRasterOperator>,
98	state: MeanRasterPixelValuesOverTimeParams,
99	}
100
101	impl InitializedPlotOperator for InitializedMeanRasterPixelValuesOverTime {
102	fn query_processor(&self) -> Result<TypedPlotQueryProcessor> {	3✔
103	let input_processor = self.raster.query_processor()?;	3✔
104	let time_position = self.state.time_position;	3✔
105	let measurement = self.raster.result_descriptor().bands[0].measurement.clone(); // TODO: adjust for multibands	3✔
106	let draw_area = self.state.area;	3✔
107
108	let processor = call_on_generic_raster_processor!(input_processor, raster => {	3✔
109	MeanRasterPixelValuesOverTimeQueryProcessor { raster, time_position, measurement, draw_area }.boxed()	3✔
110	});
111
112	Ok(TypedPlotQueryProcessor::JsonVega(processor))	3✔
113	}	3✔
114
115	fn result_descriptor(&self) -> &PlotResultDescriptor {	×
116	&self.result_descriptor	×
117	}	×
118
119	fn canonic_name(&self) -> CanonicOperatorName {	×
120	self.name.clone()	×
121	}	×
122	}
123
124	/// A query processor that calculates the `TemporalRasterMeanPlot` about its input.
125	pub struct MeanRasterPixelValuesOverTimeQueryProcessor<P: Pixel> {
126	raster: Box<dyn RasterQueryProcessor<RasterType = P>>,
127	time_position: MeanRasterPixelValuesOverTimePosition,
128	measurement: Measurement,
129	draw_area: bool,
130	}
131
132	#[async_trait]
133	impl<P: Pixel> PlotQueryProcessor for MeanRasterPixelValuesOverTimeQueryProcessor<P> {
134	type OutputFormat = PlotData;
135
136	fn plot_type(&self) -> &'static str {	×
137	MEAN_RASTER_PIXEL_VALUES_OVER_TIME_NAME	×
138	}	×
139
140	async fn plot_query<'a>(	3✔
141	&'a self,	3✔
142	query: PlotQueryRectangle,	3✔
143	ctx: &'a dyn QueryContext,	3✔
144	) -> Result<Self::OutputFormat> {	3✔
145	let means = Self::calculate_means(	3✔
146	self.raster	3✔
147	.query(	3✔
148	RasterQueryRectangle::from_qrect_and_bands(&query, BandSelection::first()),	3✔
149	ctx,	3✔
150	)	3✔
NEW 151	.await?,	×
152	self.time_position,	3✔
153	)
154	.await?;	72✔
155
156	let plot = Self::generate_plot(means, self.measurement.clone(), self.draw_area)?;	3✔
157
158	let plot_data = plot.to_vega_embeddable(false)?;	3✔
159
160	Ok(plot_data)	3✔
161	}	6✔
162	}
163
164	impl<P: Pixel> MeanRasterPixelValuesOverTimeQueryProcessor<P> {
165	async fn calculate_means(	3✔
166	mut tile_stream: BoxStream<'_, Result<RasterTile2D<P>>>,	3✔
167	position: MeanRasterPixelValuesOverTimePosition,	3✔
168	) -> Result<BTreeMap<TimeInstance, MeanCalculator>> {	3✔
169	let mut means: BTreeMap<TimeInstance, MeanCalculator> = BTreeMap::new();	3✔
170
171	while let Some(tile) = tile_stream.next().await {	44,187✔
172	let tile = tile?;	44,184✔
173
174	match tile.grid_array {	44,184✔
175	geoengine_datatypes::raster::GridOrEmpty::Grid(g) => {	20✔
176	let time = Self::time_interval_projection(tile.time, position);	20✔
177	let mean = means.entry(time).or_default();	20✔
178	mean.add(g.masked_element_deref_iterator());	20✔
179	}	20✔
180	geoengine_datatypes::raster::GridOrEmpty::Empty(_) => (),	44,164✔
181	}
182	}
183
184	Ok(means)	3✔
185	}	3✔
186
187	#[inline]
188	fn time_interval_projection(	20✔
189	time_interval: TimeInterval,	20✔
190	time_position: MeanRasterPixelValuesOverTimePosition,	20✔
191	) -> TimeInstance {	20✔
192	match time_position {	20✔
193	MeanRasterPixelValuesOverTimePosition::Start => time_interval.start(),	19✔
194	MeanRasterPixelValuesOverTimePosition::Center => TimeInstance::from_millis_unchecked(	1✔
195	average_floor(time_interval.start().inner(), time_interval.end().inner()),	1✔
196	),	1✔
197	MeanRasterPixelValuesOverTimePosition::End => time_interval.end(),	×
198	}
199	}	20✔
200
201	fn generate_plot(	3✔
202	means: BTreeMap<TimeInstance, MeanCalculator>,	3✔
203	measurement: Measurement,	3✔
204	draw_area: bool,	3✔
205	) -> Result<AreaLineChart> {	3✔
206	let mut timestamps = Vec::with_capacity(means.len());	3✔
207	let mut values = Vec::with_capacity(means.len());	3✔
208
209	for (timestamp, mean_calculator) in means {	9✔
210	timestamps.push(timestamp);	6✔
211	values.push(mean_calculator.mean());	6✔
212	}	6✔
213
214	AreaLineChart::new(timestamps, values, measurement, draw_area).map_err(Into::into)	3✔
215	}	3✔
216	}
217
218	struct MeanCalculator {
219	mean: f64,
220	n: usize,
221	}
222
223	impl Default for MeanCalculator {
224	fn default() -> Self {	6✔
225	Self { mean: 0., n: 0 }	6✔
226	}	6✔
227	}
228
229	impl MeanCalculator {
230	#[inline]
231	fn add<P: Pixel, I: Iterator<Item = Option<P>>>(&mut self, values: I) {	20✔
232	values.flatten().for_each(\|value\| {	10,024✔
233	self.add_single_value(value);	10,024✔
234	});	10,024✔
235	}	20✔
236
237	#[inline]
238	fn add_single_value<P: Pixel>(&mut self, value: P) {	10,024✔
239	let value: f64 = value.as_();	10,024✔
240		10,024✔
241	if value.is_nan() {	10,024✔
242	return;	×
243	}	10,024✔
244		10,024✔
245	self.n += 1;	10,024✔
246	let delta = value - self.mean;	10,024✔
247	self.mean += delta / (self.n as f64);	10,024✔
248	}	10,024✔
249
250	#[inline]
251	fn mean(&self) -> f64 {	6✔
252	self.mean	6✔
253	}	6✔
254	}
255
256	#[cfg(test)]
257	mod tests {
258	use super::*;
259
260	use crate::{
261	engine::{
262	ChunkByteSize, MockExecutionContext, MockQueryContext, RasterBandDescriptors,
263	RasterOperator, RasterResultDescriptor,
264	},
265	source::GdalSource,
266	};
267	use crate::{
268	mock::{MockRasterSource, MockRasterSourceParams},
269	source::GdalSourceParameters,
270	};
271	use geoengine_datatypes::primitives::{
272	BoundingBox2D, CacheHint, Measurement, PlotSeriesSelection, SpatialResolution, TimeInterval,
273	};
274	use geoengine_datatypes::{dataset::NamedData, plots::PlotMetaData, primitives::DateTime};
275	use geoengine_datatypes::{raster::TilingSpecification, spatial_reference::SpatialReference};
276	use geoengine_datatypes::{
277	raster::{Grid2D, RasterDataType, TileInformation},
278	util::test::TestDefault,
279	};
280	use serde_json::{json, Value};
281
282	#[test]	1✔
283	fn serialization() {	1✔
284	let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {	1✔
285	params: MeanRasterPixelValuesOverTimeParams {	1✔
286	time_position: MeanRasterPixelValuesOverTimePosition::Start,	1✔
287	area: true,	1✔
288	},	1✔
289	sources: SingleRasterSource {	1✔
290	raster: GdalSource {	1✔
291	params: GdalSourceParameters {	1✔
292	data: NamedData::with_system_name("test"),	1✔
293	},	1✔
294	}	1✔
295	.boxed(),	1✔
296	},	1✔
297	};	1✔
298		1✔
299	let serialized = json!({	1✔
300	"type": "MeanRasterPixelValuesOverTime",	1✔
301	"params": {	1✔
302	"timePosition": "start",	1✔
303	"area": true,	1✔
304	},	1✔
305	"sources": {	1✔
306	"raster": {	1✔
307	"type": "GdalSource",	1✔
308	"params": {	1✔
309	"data": "test"	1✔
310	}	1✔
311	}	1✔
312	},	1✔
313	"vectorSources": [],	1✔
314	})	1✔
315	.to_string();	1✔
316		1✔
317	serde_json::from_str::<Box<dyn PlotOperator>>(&serialized).unwrap();	1✔
318		1✔
319	let deserialized: MeanRasterPixelValuesOverTime =	1✔
320	serde_json::from_str(&serialized).unwrap();	1✔
321		1✔
322	assert_eq!(deserialized.params, temporal_raster_mean_plot.params);	1✔
323	}	1✔
324
325	#[tokio::test]	1✔
326	async fn single_raster() {	1✔
327	let tile_size_in_pixels = [3, 2].into();	1✔
328	let tiling_specification = TilingSpecification {	1✔
329	origin_coordinate: [0.0, 0.0].into(),	1✔
330	tile_size_in_pixels,	1✔
331	};	1✔
332	let execution_context = MockExecutionContext::new_with_tiling_spec(tiling_specification);	1✔
333		1✔
334	let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {	1✔
335	params: MeanRasterPixelValuesOverTimeParams {	1✔
336	time_position: MeanRasterPixelValuesOverTimePosition::Center,	1✔
337	area: true,	1✔
338	},	1✔
339	sources: SingleRasterSource {	1✔
340	raster: generate_mock_raster_source(	1✔
341	vec![TimeInterval::new(	1✔
342	TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),	1✔
343	TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),	1✔
344	)	1✔
345	.unwrap()],	1✔
346	vec![vec![1, 2, 3, 4, 5, 6]],	1✔
347	),	1✔
348	},	1✔
349	};	1✔
350
351	let temporal_raster_mean_plot = temporal_raster_mean_plot	1✔
352	.boxed()	1✔
353	.initialize(WorkflowOperatorPath::initialize_root(), &execution_context)	1✔
354	.await	×
355	.unwrap();	1✔
356		1✔
357	let processor = temporal_raster_mean_plot	1✔
358	.query_processor()	1✔
359	.unwrap()	1✔
360	.json_vega()	1✔
361	.unwrap();	1✔
362
363	let result = processor	1✔
364	.plot_query(	1✔
365	PlotQueryRectangle {	1✔
366	spatial_bounds: BoundingBox2D::new((-180., -90.).into(), (180., 90.).into())	1✔
367	.unwrap(),	1✔
368	time_interval: TimeInterval::default(),	1✔
369	spatial_resolution: SpatialResolution::one(),	1✔
370	attributes: PlotSeriesSelection::all(),	1✔
371	},	1✔
372	&MockQueryContext::new(ChunkByteSize::MIN),	1✔
373	)	1✔
374	.await	×
375	.unwrap();	1✔
376
377	assert!(matches!(result.metadata, PlotMetaData::None));	1✔
378
379	let vega_json: Value = serde_json::from_str(&result.vega_string).unwrap();	1✔
380		1✔
381	assert_eq!(	1✔
382	vega_json,	1✔
383	json!({	1✔
384	"$schema": "https://vega.github.io/schema/vega-lite/v4.17.0.json",	1✔
385	"data": {	1✔
386	"values": [{	1✔
387	"x": "1995-01-01T00:00:00+00:00",	1✔
388	"y": 3.5	1✔
389	}]	1✔
390	},	1✔
391	"description": "Area Plot",	1✔
392	"encoding": {	1✔
393	"x": {	1✔
394	"field": "x",	1✔
395	"title": "Time",	1✔
396	"type": "temporal"	1✔
397	},	1✔
398	"y": {	1✔
399	"field": "y",	1✔
400	"title": "",	1✔
401	"type": "quantitative"	1✔
402	}	1✔
403	},	1✔
404	"mark": {	1✔
405	"type": "area",	1✔
406	"line": true,	1✔
407	"point": true	1✔
408	}	1✔
409	})	1✔
410	);	1✔
411	}
412
413	fn generate_mock_raster_source(	2✔
414	time_intervals: Vec<TimeInterval>,	2✔
415	values_vec: Vec<Vec<u8>>,	2✔
416	) -> Box<dyn RasterOperator> {	2✔
417	assert_eq!(time_intervals.len(), values_vec.len());	2✔
418	assert!(values_vec.iter().all(\|v\| v.len() == 6));	4✔
419
420	let mut tiles = Vec::with_capacity(time_intervals.len());	2✔
421	for (time_interval, values) in time_intervals.into_iter().zip(values_vec) {	4✔
422	tiles.push(RasterTile2D::new_with_tile_info(	4✔
423	time_interval,	4✔
424	TileInformation {	4✔
425	global_geo_transform: TestDefault::test_default(),	4✔
426	global_tile_position: [0, 0].into(),	4✔
427	tile_size_in_pixels: [3, 2].into(),	4✔
428	},	4✔
429	0,	4✔
430	Grid2D::new([3, 2].into(), values).unwrap().into(),	4✔
431	CacheHint::default(),	4✔
432	));	4✔
433	}	4✔
434
435	MockRasterSource {	2✔
436	params: MockRasterSourceParams {	2✔
437	data: tiles,	2✔
438	result_descriptor: RasterResultDescriptor {	2✔
439	data_type: RasterDataType::U8,	2✔
440	spatial_reference: SpatialReference::epsg_4326().into(),	2✔
441	time: None,	2✔
442	bbox: None,	2✔
443	resolution: None,	2✔
444	bands: RasterBandDescriptors::new_single_band(),	2✔
445	},	2✔
446	},	2✔
447	}	2✔
448	.boxed()	2✔
449	}	2✔
450
451	#[tokio::test]	1✔
452	async fn raster_series() {	1✔
453	let tile_size_in_pixels = [3, 2].into();	1✔
454	let tiling_specification = TilingSpecification {	1✔
455	origin_coordinate: [0.0, 0.0].into(),	1✔
456	tile_size_in_pixels,	1✔
457	};	1✔
458	let execution_context = MockExecutionContext::new_with_tiling_spec(tiling_specification);	1✔
459		1✔
460	let temporal_raster_mean_plot = MeanRasterPixelValuesOverTime {	1✔
461	params: MeanRasterPixelValuesOverTimeParams {	1✔
462	time_position: MeanRasterPixelValuesOverTimePosition::Start,	1✔
463	area: true,	1✔
464	},	1✔
465		1✔
466	sources: SingleRasterSource {	1✔
467	raster: generate_mock_raster_source(	1✔
468	vec![	1✔
469	TimeInterval::new(	1✔
470	TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),	1✔
471	TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),	1✔
472	)	1✔
473	.unwrap(),	1✔
474	TimeInterval::new(	1✔
475	TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),	1✔
476	TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),	1✔
477	)	1✔
478	.unwrap(),	1✔
479	TimeInterval::new(	1✔
480	TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0)),	1✔
481	TimeInstance::from(DateTime::new_utc(2005, 1, 1, 0, 0, 0)),	1✔
482	)	1✔
483	.unwrap(),	1✔
484	],	1✔
485	vec![	1✔
486	vec![1, 2, 3, 4, 5, 6],	1✔
487	vec![9, 9, 8, 8, 8, 9],	1✔
488	vec![3, 4, 5, 6, 7, 8],	1✔
489	],	1✔
490	),	1✔
491	},	1✔
492	};	1✔
493
494	let temporal_raster_mean_plot = temporal_raster_mean_plot	1✔
495	.boxed()	1✔
496	.initialize(WorkflowOperatorPath::initialize_root(), &execution_context)	1✔
497	.await	×
498	.unwrap();	1✔
499		1✔
500	let processor = temporal_raster_mean_plot	1✔
501	.query_processor()	1✔
502	.unwrap()	1✔
503	.json_vega()	1✔
504	.unwrap();	1✔
505
506	let result = processor	1✔
507	.plot_query(	1✔
508	PlotQueryRectangle {	1✔
509	spatial_bounds: BoundingBox2D::new((-180., -90.).into(), (180., 90.).into())	1✔
510	.unwrap(),	1✔
511	time_interval: TimeInterval::default(),	1✔
512	spatial_resolution: SpatialResolution::one(),	1✔
513	attributes: PlotSeriesSelection::all(),	1✔
514	},	1✔
515	&MockQueryContext::new(ChunkByteSize::MIN),	1✔
516	)	1✔
517	.await	×
518	.unwrap();	1✔
519		1✔
520	assert_eq!(	1✔
521	result,	1✔
522	AreaLineChart::new(	1✔
523	vec![	1✔
524	TimeInstance::from(DateTime::new_utc(1990, 1, 1, 0, 0, 0)),	1✔
525	TimeInstance::from(DateTime::new_utc(1995, 1, 1, 0, 0, 0)),	1✔
526	TimeInstance::from(DateTime::new_utc(2000, 1, 1, 0, 0, 0))	1✔
527	],	1✔
528	vec![3.5, 8.5, 5.5],	1✔
529	Measurement::Unitless,	1✔
530	true,	1✔
531	)	1✔
532	.unwrap()	1✔
533	.to_vega_embeddable(false)	1✔
534	.unwrap()	1✔
535	);	1✔
536	}
537	}

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