rust-grib-decoder 0.1.1

Utilities to decode GRIB2 CCSDS/AEC (template 5.0=42) payloads and extract Section 7 payloads per message.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170

use std::env;
use std::error::Error;
use std::path::Path;

use rust_grib_decoder::{
    find_candidates_for_param, read_first_grib2_latlon_grid, read_grib2_section5_template42_params_by_message,
    decode_template42_try_message_payload, summarize_file,
};

#[cfg(feature = "grib-support")]
use grib::{from_reader, Grib2SubmessageDecoder};

fn main() -> Result<(), Box<dyn Error>> {
    let args: Vec<_> = env::args_os().skip(1).collect();
    if args.len() < 4 {
        eprintln!("usage: extract_param <file.grib2> <param> <lon> <lat>");
        eprintln!("example: extract_param data.grib2 msl 113.363 22.962");
        return Ok(());
    }

    let file_s = args[0].to_string_lossy().to_string();
    let path = Path::new(&file_s);
    let param = args[1].to_string_lossy().to_string();
    let lon: f64 = args[2].to_string_lossy().parse()?;
    let lat: f64 = args[3].to_string_lossy().parse()?;

    println!("Summarizing file: {}", path.display());
    let fs = summarize_file(path)?;
    for m in fs.messages.iter().take(10) {
        println!("msg {}: key='{}' repr={:?} payload_len={} decode_ok={}", m.message_index, m.key, m.repr_template, m.payload_len, m.decode_ok);
    }

    println!("\nFinding candidates for param '{}':", param);
    let candidates = find_candidates_for_param(path, &param)?;

    // Read grid early so fallback paths can use it
    let grid = read_first_grib2_latlon_grid(path)?;

    if candidates.is_empty() {
        println!("No candidates found for '{}'. Trying filename prefix fallback...", param);
        // Try filename prefix fallback using the ParamAlias file_prefixes
        let alias = rust_grib_decoder::param_alias::ParamAlias::for_param(&param);
        if let Some(fname) = path.file_name().and_then(|s| s.to_str()) {
            for prefix in alias.file_prefixes.iter() {
                if fname.starts_with(prefix) {
                    println!("Filename '{}' starts with prefix '{}', attempting message 0 decode", fname, prefix);
                    // Try to decode message 0 using template42 Section5 params if present
                    let params_vec = read_grib2_section5_template42_params_by_message(path)?;
                    if let Some(Some(params)) = params_vec.get(0) {
                        if let Ok(values) = decode_template42_try_message_payload(path, 0, grid.ni, grid.nj, params) {
                            let val = interp_value(&values, grid.ni, grid.nj, grid.lon1_deg, grid.lat1_deg, grid.di_deg, grid.dj_deg, lon, lat);
                            println!("Decoded {} samples from message 0 — interpolated value at {}x{} = {}", values.len(), lon, lat, val);
                            return Ok(());
                        } else {
                            eprintln!("template42 rust-aec decode failed for message 0");
                        }
                    }
                    #[cfg(feature = "grib-support")]
                    {
                        use grib::{from_reader, Grib2SubmessageDecoder};
                        use std::fs::File;
                        use std::io::BufReader;
                        let f = File::open(path)?;
                        let r = BufReader::new(f);
                        let grib2 = from_reader(r)?;
                        for ((midx, _sidx), sub) in grib2.iter() {
                            if midx == 0 {
                                let (w, h) = sub.grid_shape()?;
                                let decoder = Grib2SubmessageDecoder::from(sub)?;
                                match decoder.dispatch() {
                                    Ok(decoded) => {
                                        let vals: Vec<f32> = decoded.collect();
                                        let val = interp_value(&vals, w, h, grid.lon1_deg, grid.lat1_deg, grid.di_deg, grid.dj_deg, lon, lat);
                                        println!("grib crate decoded {} samples (msg 0) — interpolated value = {}", vals.len(), val);
                                        return Ok(());
                                    }
                                    Err(e) => eprintln!("grib crate decode error (msg 0): {}", e),
                                }
                            }
                        }
                    }
                }
            }
        }

        println!("No decode succeeded using filename prefix fallback.");
        return Ok(());
    }

    for (idx, reason, key, ft) in candidates.iter().take(3) {
        println!("candidate: msg={} reason={:?} key='{}' ft={:?}", idx, reason, key, ft);
    }

    // Choose the first candidate and attempt to decode it
    let (msg_idx, _, key, _ft) = &candidates[0];
    println!("\nAttempting to decode candidate msg {} ('{}')", msg_idx, key);

    let grid = read_first_grib2_latlon_grid(path)?;
    println!("Grid: {}x{} lon1={} lat1={} di={} dj={}", grid.ni, grid.nj, grid.lon1_deg, grid.lat1_deg, grid.di_deg, grid.dj_deg);

    // Try template-42 rust-aec decode using Section5 params if available
    let params_vec = read_grib2_section5_template42_params_by_message(path)?;
    if let Some(Some(params)) = params_vec.get(*msg_idx) {
        println!("Found template42 params for msg {} — trying rust-aec decode", msg_idx);
        match decode_template42_try_message_payload(path, *msg_idx, grid.ni, grid.nj, params) {
            Ok(values) => {
                let val = interp_value(&values, grid.ni, grid.nj, grid.lon1_deg, grid.lat1_deg, grid.di_deg, grid.dj_deg, lon, lat);
                println!("Decoded {} samples — interpolated value at {}x{} = {}", values.len(), lon, lat, val);
                return Ok(());
            }
            Err(e) => eprintln!("rust-aec decode failed: {}", e),
        }
    }

    #[cfg(feature = "grib-support")]
    {
        println!("Attempting to decode using `grib` crate for msg {}", msg_idx);
        use std::fs::File;
        use std::io::BufReader;
        let f = File::open(path)?;
        let r = BufReader::new(f);
        let grib2 = from_reader(r)?;
        for ((midx, _sidx), sub) in grib2.iter() {
            if midx == *msg_idx {
                let (w, h) = sub.grid_shape()?;
                let decoder = Grib2SubmessageDecoder::from(sub)?;
                match decoder.dispatch() {
                    Ok(decoded) => {
                        let vals: Vec<f32> = decoded.collect();
                        let val = interp_value(&vals, w, h, grid.lon1_deg, grid.lat1_deg, grid.di_deg, grid.dj_deg, lon, lat);
                        println!("grib crate decoded {} samples — interpolated value = {}", vals.len(), val);
                        return Ok(());
                    }
                    Err(e) => eprintln!("grib crate decode error: {}", e),
                }
            }
        }
    }

    println!("Could not decode the selected message (no template42 params and grib-support not enabled or decoding failed).");
    Ok(())
}

fn interp_value(values: &[f32], width: usize, height: usize, lon1: f64, lat1: f64, di: f64, dj: f64, lon: f64, lat: f64) -> f64 {
    let mut dlon = lon - lon1;
    while dlon < 0.0 { dlon += 360.0; }
    while dlon >= 360.0 { dlon -= 360.0; }
    let i_f = dlon / di;
    let j_f = (lat1 - lat) / dj;
    if !i_f.is_finite() || !j_f.is_finite() { return f64::NAN; }
    let ni = width as isize; let nj = height as isize;

    let i0 = i_f.floor() as isize; let j0 = j_f.floor() as isize; let i1 = i0 + 1; let j1 = j0 + 1;
    let dx = i_f - (i0 as f64); let dy = j_f - (j0 as f64);
    let i0w = (((i0 % ni) + ni) % ni) as usize; let i1w = (((i1 % ni) + ni) % ni) as usize;
    let j0c = j0.clamp(0, nj - 1) as usize; let j1c = j1.clamp(0, nj - 1) as usize;
    let idx00 = j0c.saturating_mul(width).saturating_add(i0w);
    let idx10 = j0c.saturating_mul(width).saturating_add(i1w);
    let idx01 = j1c.saturating_mul(width).saturating_add(i0w);
    let idx11 = j1c.saturating_mul(width).saturating_add(i1w);
    let s00 = values.get(idx00).copied().unwrap_or(f32::NAN);
    let s10 = values.get(idx10).copied().unwrap_or(f32::NAN);
    let s01 = values.get(idx01).copied().unwrap_or(f32::NAN);
    let s11 = values.get(idx11).copied().unwrap_or(f32::NAN);
    if !s00.is_finite() || !s10.is_finite() || !s01.is_finite() || !s11.is_finite() { return f64::NAN; }
    let v0 = s00 * (1.0 - dx as f32) + s10 * (dx as f32);
    let v1 = s01 * (1.0 - dx as f32) + s11 * (dx as f32);
    let v = v0 * (1.0 - dy as f32) + v1 * (dy as f32);
    v as f64
}