Skip to main content

grib/
grid.rs

1#[cfg(feature = "gridpoints-proj")]
2use helpers::ProjectionLatLonIterator;
3use helpers::RegularGridIterator;
4
5pub use self::{gaussian::compute_gaussian_latitudes, rotation::Unrotate};
6use crate::{
7    GribError, GridDefinition, TryFromSlice,
8    def::grib2::template::{
9        Template3_0, Template3_1, Template3_10, Template3_20, Template3_30, Template3_40,
10        Template3_41,
11        param_set::{Grid, ScanningMode},
12    },
13};
14
15#[derive(Debug, PartialEq)]
16#[non_exhaustive]
17pub enum GridDefinitionTemplateValues {
18    Template0(Template3_0),
19    Template1(Template3_1),
20    Template10(Template3_10),
21    Template20(Template3_20),
22    Template30(Template3_30),
23    Template40(Template3_40),
24    Template41(Template3_41),
25}
26
27impl GridShortName for GridDefinitionTemplateValues {
28    fn short_name(&self) -> &'static str {
29        match self {
30            Self::Template0(def) => def.lat_lon.short_name(),
31            Self::Template1(def) => def.short_name(),
32            Self::Template10(def) => def.short_name(),
33            Self::Template20(def) => def.short_name(),
34            Self::Template30(def) => def.short_name(),
35            Self::Template40(def) => def.gaussian.short_name(),
36            Self::Template41(def) => def.short_name(),
37        }
38    }
39}
40
41impl GridPointIndex for GridDefinitionTemplateValues {
42    fn grid_shape(&self) -> (usize, usize) {
43        match self {
44            Self::Template0(def) => def.lat_lon.grid_shape(),
45            Self::Template1(def) => def.grid_shape(),
46            Self::Template10(def) => def.grid_shape(),
47            Self::Template20(def) => def.grid_shape(),
48            Self::Template30(def) => def.grid_shape(),
49            Self::Template40(def) => def.gaussian.grid_shape(),
50            Self::Template41(def) => def.grid_shape(),
51        }
52    }
53
54    fn scanning_mode(&self) -> &ScanningMode {
55        match self {
56            Self::Template0(def) => def.lat_lon.scanning_mode(),
57            Self::Template1(def) => def.scanning_mode(),
58            Self::Template10(def) => def.scanning_mode(),
59            Self::Template20(def) => def.scanning_mode(),
60            Self::Template30(def) => def.scanning_mode(),
61            Self::Template40(def) => def.gaussian.scanning_mode(),
62            Self::Template41(def) => def.scanning_mode(),
63        }
64    }
65}
66
67impl LatLons for GridDefinitionTemplateValues {
68    type Iter<'a>
69        = GridPointLatLons
70    where
71        Self: 'a;
72
73    fn latlons_unchecked<'a>(&'a self) -> Result<Self::Iter<'a>, GribError> {
74        let iter = match self {
75            Self::Template0(def) => GridPointLatLons::from(def.lat_lon.latlons_unchecked()?),
76            Self::Template1(def) => GridPointLatLons::from(def.latlons_unchecked()?),
77            Self::Template10(def) => GridPointLatLons::from(def.latlons_unchecked()?),
78            #[cfg(feature = "gridpoints-proj")]
79            Self::Template20(def) => GridPointLatLons::from(def.latlons_unchecked()?),
80            Self::Template30(def) => GridPointLatLons::from(def.latlons_unchecked()?),
81            Self::Template40(def) => GridPointLatLons::from(def.gaussian.latlons_unchecked()?),
82            Self::Template41(def) => GridPointLatLons::from(def.latlons_unchecked()?),
83            #[cfg(not(feature = "gridpoints-proj"))]
84            _ => {
85                return Err(GribError::NotSupported(
86                    "lat/lon computation support for the template is dropped in this build"
87                        .to_owned(),
88                ));
89            }
90        };
91        Ok(iter)
92    }
93}
94
95impl TryFrom<&GridDefinition> for GridDefinitionTemplateValues {
96    type Error = GribError;
97
98    fn try_from(value: &GridDefinition) -> Result<Self, Self::Error> {
99        // In the future, we should switch the implementation like this:
100        //
101        // ```
102        // let buf = &value.payload;
103        // let mut pos = 0;
104        // let payload = crate::def::grib2::Section3Payload::try_from_slice(buf, &mut pos)?;
105        // let template = match payload.template {
106        // ..
107        // }
108        // ```
109        //
110        // However, since the current implementation of the templates has many
111        // limitations, to prevent errors, the template reading process is
112        // implemented as follows.
113        let buf = &value.payload[9..];
114        let mut pos = 0;
115        let num = value.grid_tmpl_num();
116        let template = match num {
117            0 => {
118                GridDefinitionTemplateValues::Template0(Template3_0::try_from_slice(buf, &mut pos)?)
119            }
120            1 => {
121                GridDefinitionTemplateValues::Template1(Template3_1::try_from_slice(buf, &mut pos)?)
122            }
123            10 => GridDefinitionTemplateValues::Template10(Template3_10::try_from_slice(
124                buf, &mut pos,
125            )?),
126            20 => GridDefinitionTemplateValues::Template20(Template3_20::try_from_slice(
127                buf, &mut pos,
128            )?),
129            30 => GridDefinitionTemplateValues::Template30(Template3_30::try_from_slice(
130                buf, &mut pos,
131            )?),
132            40 => GridDefinitionTemplateValues::Template40(Template3_40::try_from_slice(
133                buf, &mut pos,
134            )?),
135            41 => GridDefinitionTemplateValues::Template41(Template3_41::try_from_slice(
136                buf, &mut pos,
137            )?),
138            _ => {
139                return Err(GribError::NotSupported(format!(
140                    "lat/lon computation support for the template {num} is dropped in this build"
141                )));
142            }
143        };
144        if buf.len() > pos {
145            return Err(GribError::NotSupported(
146                "template with list of number of points".to_owned(),
147            ));
148        }
149        Ok(template)
150    }
151}
152
153/// A functionality to return a short name of the grid system.
154pub trait GridShortName {
155    /// Returns the grid type.
156    ///
157    /// The grid types are denoted as short strings based on `gridType` used in
158    /// ecCodes.
159    ///
160    /// This is provided primarily for debugging and simple notation purposes.
161    /// It is better to use enum variants instead of the string notation to
162    /// determine the grid type.
163    fn short_name(&self) -> &'static str;
164}
165
166impl<T: GridShortName + ?Sized> GridShortName for &T {
167    fn short_name(&self) -> &'static str {
168        (*self).short_name()
169    }
170}
171
172/// A functionality to generate an iterator over latitude/longitude of grid
173/// points.
174pub trait LatLons {
175    type Iter<'a>: Iterator<Item = (f32, f32)>
176    where
177        Self: 'a;
178
179    /// Computes and returns an iterator over latitudes and longitudes of grid
180    /// points in degrees. Unlike the return values from [`LatLons::latlons`],
181    /// the longitude values from thie method do not necessarily fall within
182    /// the range `[-180°, 180°]`.
183    fn latlons_unchecked<'a>(&'a self) -> Result<Self::Iter<'a>, GribError>;
184
185    /// Computes and returns an iterator over latitudes and longitudes of grid
186    /// points in degrees. The returned longitude values are converted to fall
187    /// within the range `[-180°, 180°]`.
188    ///
189    /// The order of lat/lon data of grid points is the same as the order of the
190    /// grid point values, defined by the scanning mode
191    /// ([`ScanningMode`](`crate::def::grib2::template::param_set::ScanningMode`)) in the data.
192    #[allow(clippy::type_complexity)]
193    fn latlons<'a>(
194        &'a self,
195    ) -> Result<std::iter::Map<Self::Iter<'a>, fn((f32, f32)) -> (f32, f32)>, GribError> {
196        let iter = self
197            .latlons_unchecked()?
198            .map(helpers::normalize_latlon as fn((f32, f32)) -> (f32, f32));
199        Ok(iter)
200    }
201}
202
203impl<T: LatLons + ?Sized> LatLons for &T {
204    type Iter<'a>
205        = <T as LatLons>::Iter<'a>
206    where
207        Self: 'a;
208
209    fn latlons_unchecked<'a>(&'a self) -> Result<Self::Iter<'a>, GribError> {
210        (*self).latlons_unchecked()
211    }
212}
213
214/// An iterator over latitudes and longitudes of grid points in a submessage.
215///
216/// This `struct` is created by the [`latlons`] method on [`LatLons`]
217/// implemented for [`SubMessage`]. See its documentation for more.
218///
219/// [`latlons`]: crate::context::SubMessage::latlons
220/// [`SubMessage`]: crate::context::SubMessage
221#[derive(Clone)]
222pub struct GridPointLatLons(LatLonsWrapper);
223
224impl Iterator for GridPointLatLons {
225    type Item = (f32, f32);
226
227    fn next(&mut self) -> Option<Self::Item> {
228        match self {
229            Self(LatLonsWrapper::SigR(iter)) => iter.next(),
230            Self(LatLonsWrapper::SigUR(iter)) => iter.next(),
231            Self(LatLonsWrapper::SigIf(iter)) => iter.next(),
232            #[cfg(feature = "gridpoints-proj")]
233            Self(LatLonsWrapper::SigP(iter)) => iter.next(),
234        }
235    }
236
237    fn size_hint(&self) -> (usize, Option<usize>) {
238        match self {
239            Self(LatLonsWrapper::SigR(iter)) => iter.size_hint(),
240            Self(LatLonsWrapper::SigUR(iter)) => iter.size_hint(),
241            Self(LatLonsWrapper::SigIf(iter)) => iter.size_hint(),
242            #[cfg(feature = "gridpoints-proj")]
243            Self(LatLonsWrapper::SigP(iter)) => iter.size_hint(),
244        }
245    }
246}
247
248impl From<RegularGridIterator> for GridPointLatLons {
249    fn from(value: RegularGridIterator) -> Self {
250        Self(LatLonsWrapper::SigR(value))
251    }
252}
253
254impl From<Unrotate<RegularGridIterator>> for GridPointLatLons {
255    fn from(value: Unrotate<RegularGridIterator>) -> Self {
256        Self(LatLonsWrapper::SigUR(value))
257    }
258}
259
260impl From<std::vec::IntoIter<(f32, f32)>> for GridPointLatLons {
261    fn from(value: std::vec::IntoIter<(f32, f32)>) -> Self {
262        Self(LatLonsWrapper::SigIf(value))
263    }
264}
265
266#[cfg(feature = "gridpoints-proj")]
267impl From<ProjectionLatLonIterator> for GridPointLatLons {
268    fn from(value: ProjectionLatLonIterator) -> Self {
269        Self(LatLonsWrapper::SigP(value))
270    }
271}
272
273#[derive(Clone)]
274enum LatLonsWrapper {
275    SigR(RegularGridIterator),
276    SigUR(Unrotate<RegularGridIterator>),
277    SigIf(std::vec::IntoIter<(f32, f32)>),
278    #[cfg(feature = "gridpoints-proj")]
279    SigP(ProjectionLatLonIterator),
280}
281
282/// A functionality to generate an iterator over 2D index `(i, j)` of grid
283/// points.
284///
285/// # Examples
286///
287/// ```
288/// use grib::{GridPointIndex, def::grib2::template::param_set::ScanningMode};
289///
290/// struct Grid {
291///     ni: u32,
292///     nj: u32,
293///     scanning_mode: ScanningMode,
294/// }
295///
296/// impl GridPointIndex for Grid {
297///     fn grid_shape(&self) -> (usize, usize) {
298///         (self.ni as usize, self.nj as usize)
299///     }
300///
301///     fn scanning_mode(&self) -> &ScanningMode {
302///         &self.scanning_mode
303///     }
304/// }
305///
306/// let grid_2x3 = Grid {
307///     ni: 2,
308///     nj: 3,
309///     scanning_mode: ScanningMode(0b01000000),
310/// };
311/// assert_eq!(grid_2x3.grid_shape(), (2, 3));
312///
313/// let mut iter = grid_2x3.ij().unwrap();
314/// assert_eq!(iter.next(), Some((0, 0)));
315/// assert_eq!(iter.next(), Some((1, 0)));
316/// assert_eq!(iter.next(), Some((0, 1)));
317/// assert_eq!(iter.next(), Some((1, 1)));
318/// assert_eq!(iter.next(), Some((0, 2)));
319/// assert_eq!(iter.next(), Some((1, 2)));
320/// assert_eq!(iter.next(), None);
321/// ```
322pub trait GridPointIndex {
323    /// Returns the shape of the grid, i.e. a tuple of the number of grids in
324    /// the i and j directions.
325    fn grid_shape(&self) -> (usize, usize);
326
327    /// Returns [`ScanningMode`] used for the iteration.
328    fn scanning_mode(&self) -> &ScanningMode;
329
330    /// Returns an iterator over 2D index `(i, j)` of grid points.
331    fn ij(&self) -> Result<GridPointIndexIterator, GribError> {
332        GridPointIndexIterator::new(self.grid_shape(), *self.scanning_mode())
333    }
334}
335
336impl<T: GridPointIndex + ?Sized> GridPointIndex for &T {
337    fn grid_shape(&self) -> (usize, usize) {
338        (*self).grid_shape()
339    }
340
341    fn scanning_mode(&self) -> &ScanningMode {
342        (*self).scanning_mode()
343    }
344}
345
346/// An iterator over 2D index `(i, j)` of grid points.
347///
348/// This `struct` is created by the [`GridPointIndex::ij`] method. See its
349/// documentation for more.
350#[derive(Clone)]
351pub struct GridPointIndexIterator {
352    major_len: usize,
353    minor_len: usize,
354    scanning_mode: ScanningMode,
355    major_pos: usize,
356    minor_pos: usize,
357    increments: bool,
358}
359
360impl GridPointIndexIterator {
361    pub(crate) fn new(
362        (i_len, j_len): (usize, usize),
363        scanning_mode: ScanningMode,
364    ) -> Result<Self, GribError> {
365        if scanning_mode.has_unsupported_flags() {
366            let ScanningMode(mode) = scanning_mode;
367            return Err(GribError::NotSupported(format!("scanning mode {mode}")));
368        }
369
370        let (major_len, minor_len) = if scanning_mode.is_consecutive_for_i() {
371            (j_len, i_len)
372        } else {
373            (i_len, j_len)
374        };
375
376        Ok(Self {
377            major_len,
378            minor_len,
379            scanning_mode,
380            minor_pos: 0,
381            major_pos: 0,
382            increments: true,
383        })
384    }
385}
386
387impl Iterator for GridPointIndexIterator {
388    type Item = (usize, usize);
389
390    fn next(&mut self) -> Option<Self::Item> {
391        if self.major_pos == self.major_len {
392            return None;
393        }
394
395        let minor = if self.increments {
396            self.minor_pos
397        } else {
398            self.minor_len - self.minor_pos - 1
399        };
400        let major = self.major_pos;
401
402        self.minor_pos += 1;
403        if self.minor_pos == self.minor_len {
404            self.major_pos += 1;
405            self.minor_pos = 0;
406            if self.scanning_mode.scans_alternating_rows() {
407                self.increments = !self.increments;
408            }
409        }
410
411        if self.scanning_mode.is_consecutive_for_i() {
412            Some((minor, major))
413        } else {
414            Some((major, minor))
415        }
416    }
417
418    fn size_hint(&self) -> (usize, Option<usize>) {
419        let len = (self.major_len - self.major_pos) * self.minor_len - self.minor_pos;
420        (len, Some(len))
421    }
422}
423
424pub(crate) trait AngleUnit {
425    fn angle_unit(&self) -> f64;
426}
427
428impl AngleUnit for Grid {
429    fn angle_unit(&self) -> f64 {
430        let basic_angle = self.initial_production_domain_basic_angle;
431        let sub_angle = self.basic_angle_subdivisions;
432        if basic_angle == 0 {
433            1e-6
434        } else {
435            basic_angle as f64 / sub_angle as f64
436        }
437    }
438}
439
440#[cfg(test)]
441mod tests {
442    use super::*;
443    use crate::def::grib2::template::param_set::{EarthShape, ScaledValue};
444
445    #[test]
446    fn grid_definition_template_0() -> Result<(), Box<dyn std::error::Error>> {
447        let buf = crate::test_utils::decompress_to_vec(
448            crate::test_utils::data::grib2::JMA_TORNADO_NOWCAST,
449        )?;
450        let data =
451            GridDefinition::from_payload(buf[0x2a..0x6d].to_vec().into_boxed_slice()).unwrap();
452
453        let actual = GridDefinitionTemplateValues::try_from(&data).unwrap();
454        let expected = GridDefinitionTemplateValues::Template0(Template3_0 {
455            earth: EarthShape {
456                shape: 4,
457                spherical_earth_radius: ScaledValue {
458                    scale_factor: 0xff,
459                    scaled_value: 0xffffffff,
460                },
461                major_axis: ScaledValue {
462                    scale_factor: 1,
463                    scaled_value: 63781370,
464                },
465                minor_axis: ScaledValue {
466                    scale_factor: 1,
467                    scaled_value: 63567523,
468                },
469            },
470            lat_lon: crate::def::grib2::template::param_set::LatLonGrid {
471                grid: crate::def::grib2::template::param_set::Grid {
472                    ni: 256,
473                    nj: 336,
474                    initial_production_domain_basic_angle: 0,
475                    basic_angle_subdivisions: 0xffffffff,
476                    first_point_lat: 47958333,
477                    first_point_lon: 118062500,
478                    resolution_and_component_flags:
479                        crate::def::grib2::template::param_set::ResolutionAndComponentFlags(
480                            0b00110000,
481                        ),
482                    last_point_lat: 20041667,
483                    last_point_lon: 149937500,
484                },
485                i_direction_inc: 125000,
486                j_direction_inc: 83333,
487                scanning_mode: crate::def::grib2::template::param_set::ScanningMode(0b00000000),
488            },
489        });
490        assert_eq!(actual, expected);
491        Ok(())
492    }
493}
494
495mod earth;
496mod flags;
497mod gaussian;
498mod helpers;
499mod lambert;
500mod latlon;
501mod mercator;
502mod polar_stereographic;
503mod rotation;