wsi_dicom_annotations/annotations/parametric_map/
preview.rs

1use crate::{DicomCode, Error, Point2, Result};
2
3use super::source::PixelTile;
4use super::ParametricMapDocument;
5
6const PREVIEW_CHUNK_SIDE: u32 = 256;
7const MAX_PREVIEW_PIXELS: usize = 1_048_576;
8
9#[derive(Debug, Clone)]
10pub struct ParametricMapPreview {
11    width: u32,
12    height: u32,
13    normalized_values: Vec<f32>,
14    minimum: f64,
15    maximum: f64,
16    channel_name: String,
17    quantity: DicomCode,
18    unit: DicomCode,
19    base_corners: [Point2; 4],
20}
21
22impl ParametricMapPreview {
23    #[must_use]
24    pub const fn dimensions(&self) -> (u32, u32) {
25        (self.width, self.height)
26    }
27
28    #[must_use]
29    pub fn normalized_values(&self) -> &[f32] {
30        &self.normalized_values
31    }
32
33    #[must_use]
34    pub const fn value_range(&self) -> (f64, f64) {
35        (self.minimum, self.maximum)
36    }
37
38    #[must_use]
39    pub fn channel_name(&self) -> &str {
40        &self.channel_name
41    }
42
43    #[must_use]
44    pub const fn quantity(&self) -> &DicomCode {
45        &self.quantity
46    }
47
48    #[must_use]
49    pub const fn unit(&self) -> &DicomCode {
50        &self.unit
51    }
52
53    #[must_use]
54    pub const fn base_corners(&self) -> &[Point2; 4] {
55        &self.base_corners
56    }
57}
58
59impl ParametricMapDocument {
60    /// Build a bounded, display-only heatmap for one selected parameter.
61    ///
62    /// Values are averaged after normalization into the parameter's finite
63    /// Real World Value range. Missing bins remain NaN. This scans the source
64    /// in bounded chunks and never allocates the full raster.
65    pub fn preview(
66        &self,
67        selected_channel_ordinal: usize,
68        maximum_pixels: usize,
69    ) -> Result<ParametricMapPreview> {
70        if maximum_pixels == 0 || maximum_pixels > MAX_PREVIEW_PIXELS {
71            return Err(Error::InvalidInput(format!(
72                "PM preview pixel limit must be 1..={MAX_PREVIEW_PIXELS}"
73            )));
74        }
75        let channel = *self
76            .selected_channels
77            .get(selected_channel_ordinal)
78            .ok_or_else(|| {
79                Error::InvalidInput(format!(
80                    "selected PM preview channel ordinal {selected_channel_ordinal} is outside 0..{}",
81                    self.selected_channels.len()
82                ))
83            })?;
84        let range = self.channel_ranges[selected_channel_ordinal];
85        let (width, height) = preview_dimensions(
86            self.descriptor.width,
87            self.descriptor.height,
88            maximum_pixels,
89        )?;
90        let normalized_values = read_preview_values(self, channel, width, height, range)?;
91        let channel_profile = &self.profile.channels[channel];
92        Ok(ParametricMapPreview {
93            width,
94            height,
95            normalized_values,
96            minimum: range.minimum,
97            maximum: range.maximum,
98            channel_name: channel_profile.name.clone(),
99            quantity: channel_profile.quantity.clone(),
100            unit: channel_profile.unit.clone(),
101            base_corners: base_corners(self)?,
102        })
103    }
104}
105
106fn read_preview_values(
107    document: &ParametricMapDocument,
108    channel: usize,
109    width: u32,
110    height: u32,
111    range: super::document::ValueRange,
112) -> Result<Vec<f32>> {
113    let mut accumulator = PreviewAccumulator::new(
114        document.descriptor.width,
115        document.descriptor.height,
116        width,
117        height,
118        range,
119    )?;
120    for origin_y in (0..document.descriptor.height).step_by(PREVIEW_CHUNK_SIDE as usize) {
121        let chunk_height = PREVIEW_CHUNK_SIDE.min(document.descriptor.height - origin_y);
122        for origin_x in (0..document.descriptor.width).step_by(PREVIEW_CHUNK_SIDE as usize) {
123            let chunk_width = PREVIEW_CHUNK_SIDE.min(document.descriptor.width - origin_x);
124            let tile = document.raster.read_tile(
125                channel,
126                origin_y,
127                origin_x,
128                chunk_height,
129                chunk_width,
130            )?;
131            accumulator.add_tile(&tile, origin_y, origin_x, chunk_width);
132        }
133    }
134    Ok(accumulator.finish())
135}
136
137fn preview_dimensions(width: u32, height: u32, maximum_pixels: usize) -> Result<(u32, u32)> {
138    let source_pixels = u64::from(width) * u64::from(height);
139    if source_pixels <= maximum_pixels as u64 {
140        return Ok((width, height));
141    }
142    let aspect = f64::from(width) / f64::from(height);
143    let mut preview_width = ((maximum_pixels as f64 * aspect).sqrt().round() as u64)
144        .clamp(1, u64::from(width).min(maximum_pixels as u64));
145    let preview_height = (maximum_pixels as u64 / preview_width).clamp(1, u64::from(height));
146    preview_width = preview_width.min(maximum_pixels as u64 / preview_height);
147    if preview_width == 0 || preview_height == 0 {
148        return Err(Error::InvalidInput(
149            "PM preview dimensions resolved to zero".into(),
150        ));
151    }
152    Ok((preview_width as u32, preview_height as u32))
153}
154
155struct PreviewAccumulator {
156    source_width: u32,
157    source_height: u32,
158    preview_width: u32,
159    preview_height: u32,
160    range: super::document::ValueRange,
161    sums: Vec<f64>,
162    counts: Vec<u64>,
163}
164
165impl PreviewAccumulator {
166    fn new(
167        source_width: u32,
168        source_height: u32,
169        preview_width: u32,
170        preview_height: u32,
171        range: super::document::ValueRange,
172    ) -> Result<Self> {
173        let length = usize::try_from(u64::from(preview_width) * u64::from(preview_height))
174            .map_err(|_| {
175                Error::InvalidInput("PM preview dimensions exceed addressable memory".into())
176            })?;
177        Ok(Self {
178            source_width,
179            source_height,
180            preview_width,
181            preview_height,
182            range,
183            sums: vec![0.0; length],
184            counts: vec![0; length],
185        })
186    }
187
188    fn add_tile(&mut self, tile: &PixelTile, origin_y: u32, origin_x: u32, chunk_width: u32) {
189        match tile {
190            PixelTile::F32(values) => values.iter().enumerate().for_each(|(index, value)| {
191                self.observe(index, f64::from(*value), origin_y, origin_x, chunk_width);
192            }),
193            PixelTile::F64(values) => values.iter().enumerate().for_each(|(index, value)| {
194                self.observe(index, *value, origin_y, origin_x, chunk_width);
195            }),
196        }
197    }
198
199    fn observe(
200        &mut self,
201        index: usize,
202        value: f64,
203        origin_y: u32,
204        origin_x: u32,
205        chunk_width: u32,
206    ) {
207        if !value.is_finite() {
208            return;
209        }
210        let row = origin_y + index as u32 / chunk_width;
211        let column = origin_x + index as u32 % chunk_width;
212        let preview_row =
213            u64::from(row) * u64::from(self.preview_height) / u64::from(self.source_height);
214        let preview_column =
215            u64::from(column) * u64::from(self.preview_width) / u64::from(self.source_width);
216        let preview_index = (preview_row * u64::from(self.preview_width) + preview_column) as usize;
217        let normalized = if self.range.minimum == self.range.maximum {
218            0.5
219        } else {
220            ((value - self.range.minimum) / (self.range.maximum - self.range.minimum))
221                .clamp(0.0, 1.0)
222        };
223        self.sums[preview_index] += normalized;
224        self.counts[preview_index] += 1;
225    }
226
227    fn finish(self) -> Vec<f32> {
228        self.sums
229            .into_iter()
230            .zip(self.counts)
231            .map(|(sum, count)| {
232                if count == 0 {
233                    f32::NAN
234                } else {
235                    (sum / count as f64) as f32
236                }
237            })
238            .collect()
239    }
240}
241
242fn base_corners(document: &ParametricMapDocument) -> Result<[Point2; 4]> {
243    let right = f64::from(document.descriptor.width) - 0.5;
244    let bottom = f64::from(document.descriptor.height) - 0.5;
245    let point = |column, row| {
246        let position = document.geometry.position_at(column, row);
247        document
248            .source
249            .slide_coordinate_to_pixel3(position[0], position[1], position[2])
250    };
251    Ok([
252        point(-0.5, -0.5)?,
253        point(right, -0.5)?,
254        point(right, bottom)?,
255        point(-0.5, bottom)?,
256    ])
257}
wsi_dicom_annotations/annotations/parametric_map/preview.rs

wsi_dicom_annotations/annotations/parametric_map/
preview.rs
