diff --git a/Cargo.lock b/Cargo.lock index 45decf98..8e747a40 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -3394,6 +3394,7 @@ dependencies = [ "photocraft-ops", "photocraft-paint", "photocraft-raster", + "photocraft-raw", "photocraft-text", "photocraft-vector", "rayon", @@ -3457,6 +3458,7 @@ dependencies = [ "photocraft-geom", "photocraft-psd", "photocraft-raster", + "photocraft-raw", "photocraft-text", "photocraft-vector", "proptest", @@ -3502,6 +3504,15 @@ dependencies = [ "serde", ] +[[package]] +name = "photocraft-raw" +version = "0.1.1" +dependencies = [ + "photocraft-raw", + "png 0.18.1", + "rayon", +] + [[package]] name = "photocraft-testkit" version = "0.1.1" diff --git a/Cargo.toml b/Cargo.toml index 02d9e5c3..448a1743 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -19,6 +19,7 @@ photocraft-doc = { path = "crates/doc" } photocraft-ops = { path = "crates/ops" } photocraft-psd = { path = "crates/psd" } photocraft-codecs = { path = "crates/codecs" } +photocraft-raw = { path = "crates/raw" } photocraft-compose = { path = "crates/compose" } photocraft-gpu = { path = "crates/gpu" } photocraft-paint = { path = "crates/paint" } diff --git a/apps/photocraft-web/src/web.rs b/apps/photocraft-web/src/web.rs index 6472bacf..a010403a 100644 --- a/apps/photocraft-web/src/web.rs +++ b/apps/photocraft-web/src/web.rs @@ -12,8 +12,10 @@ use wasm_bindgen::JsCast as _; type Inbox = Arc)>>>; /// Everything File › Open reads: PhotoCraft documents, Photoshop documents and flat images. -const OPEN_EXTS: &[&str] = - &["pcraft", "psd", "psb", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "ico", "qoi", "exr", "hdr", "pbm", "pgm", "ppm", "pam", "pfm"]; +const OPEN_EXTS: &[&str] = &[ + "pcraft", "psd", "psb", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "ico", "qoi", "exr", "hdr", "pbm", "pgm", "ppm", "pam", "pfm", + "dng", "cr2", "cr3", "nef", "nrw", "arw", "pef", "orf", "rw2", "raf", +]; const CANVAS_ID: &str = "photocraft_canvas"; pub fn start() { diff --git a/apps/photocraft/src/services.rs b/apps/photocraft/src/services.rs index ef155fd5..bdbf4c30 100644 --- a/apps/photocraft/src/services.rs +++ b/apps/photocraft/src/services.rs @@ -13,8 +13,10 @@ use std::rc::Rc; use std::sync::Arc; /// Everything File › Open reads: PhotoCraft documents, Photoshop documents and flat images. -const OPEN_EXTS: &[&str] = - &["pcraft", "psd", "psb", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "ico", "qoi", "exr", "hdr", "pbm", "pgm", "ppm", "pam", "pfm"]; +const OPEN_EXTS: &[&str] = &[ + "pcraft", "psd", "psb", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "ico", "qoi", "exr", "hdr", "pbm", "pgm", "ppm", "pam", "pfm", + "dng", "cr2", "cr3", "nef", "nrw", "arw", "pef", "orf", "rw2", "raf", +]; /// File › Save As formats: (filter name, extensions). The filter matching the suggested name's /// extension comes first, so a .pcraft document saves as .pcraft by default and everything else diff --git a/crates/codecs/README.md b/crates/codecs/README.md index 88af8a54..a88bcf30 100644 --- a/crates/codecs/README.md +++ b/crates/codecs/README.md @@ -60,6 +60,10 @@ the same. ## Documented asymmetries and limitations +* **Camera raw files** (DNG, CR2, NEF, ARW… which are TIFF-structured) are recognised and refused + with `CodecError::Unsupported`: they are sensor data, not flat images. `photocraft-raw` decodes + and develops them, and `photocraft-io` routes them there. + * **AVIF (the only asymmetric format).** Encoding uses `ravif`, which is pure Rust. Decoding needs `dav1d`, which is C. AVIF is therefore read-unsupported, and write support is gated behind the non-default `avif` feature. In a default build it is neither readable nor writable, diff --git a/crates/codecs/src/codecs/tiff.rs b/crates/codecs/src/codecs/tiff.rs index ce15c2ce..03bfd323 100644 --- a/crates/codecs/src/codecs/tiff.rs +++ b/crates/codecs/src/codecs/tiff.rs @@ -135,7 +135,7 @@ fn camera_raw(b: &[u8]) -> bool { pub(crate) fn decode(bytes: &[u8], limits: &Limits) -> Result { if camera_raw(bytes) { - return Err(CodecError::unsupported(F, "camera raw files (such as CR2, NEF, ARW or DNG) are not supported")); + return Err(CodecError::unsupported(F, "camera raw files (such as CR2, NEF, ARW or DNG) are not flat images; decode them with photocraft-raw")); } let tl = { let mut l = tiff::decoder::Limits::default(); diff --git a/crates/engine/Cargo.toml b/crates/engine/Cargo.toml index 5c05f310..9f4bb571 100644 --- a/crates/engine/Cargo.toml +++ b/crates/engine/Cargo.toml @@ -28,3 +28,7 @@ photocraft-format = { workspace = true } blake3 = { version = "1.8", default-features = false, features = ["pure"] } photocraft-codecs = { workspace = true } flate2 = { workspace = true } + +[dev-dependencies] +# Synthetic camera raw files for the open tests. +photocraft-raw = { workspace = true, features = ["testgen"] } diff --git a/crates/engine/src/file_cmds.rs b/crates/engine/src/file_cmds.rs index 52ac78b9..6ec445b2 100644 --- a/crates/engine/src/file_cmds.rs +++ b/crates/engine/src/file_cmds.rs @@ -114,8 +114,10 @@ pub(crate) fn list_images(dir: &str) -> Result> { } /// Extensions the batch commands pick up from a folder. -const OPENABLE: &[&str] = - &["psd", "psb", "pcraft", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "exr", "hdr", "qoi", "ico", "pnm", "ppm", "pgm"]; +const OPENABLE: &[&str] = &[ + "psd", "psb", "pcraft", "png", "jpg", "jpeg", "tif", "tiff", "webp", "gif", "bmp", "tga", "exr", "hdr", "qoi", "ico", "pnm", "ppm", "pgm", "dng", "cr2", + "nef", "nrw", "arw", "pef", +]; pub(crate) fn file_name(path: &str) -> String { path.rsplit(['/', '\\']).next().unwrap_or(path).to_string() @@ -265,11 +267,13 @@ pub fn open_bytes_as(s: &mut Session, name: &str, bytes: &[u8], as_ext: Option<& Some(ext) => format!("{}.{}", stem(name), ext.trim_start_matches('.')), None => name.to_string(), }; - let mut doc = import(&decode_name, bytes)?; + let r = photocraft_io::import(&decode_name, bytes).map_err(|e| EngineError::Other(format!("{decode_name}: {e}")))?; + let mut doc = r.document; doc.name = file_name(name); // Color Settings policies (preserve / convert / discard the embedded profile). let (i, color) = s.open_document(doc, path); - Ok(json!({"document": i, "color": color})) + // Import notes (e.g. how a camera raw was developed, or that only its preview opened). + Ok(json!({"document": i, "color": color, "warnings": r.warnings})) } fn open_as(s: &mut Session, p: &Value) -> Result { diff --git a/crates/engine/src/file_cmds/tests.rs b/crates/engine/src/file_cmds/tests.rs index 8b69c6f3..1fa71ef6 100644 --- a/crates/engine/src/file_cmds/tests.rs +++ b/crates/engine/src/file_cmds/tests.rs @@ -48,6 +48,23 @@ fn close_all_and_close_others() { assert!(!s.is_enabled("file.closeAll")); } +#[test] +fn camera_raw_opens_as_16_bit_with_notes() { + use photocraft_raw::testgen::{DngSpec, mosaic, scene}; + let (w, h) = (24, 16); + let mut spec = DngSpec::cfa(w, h, mosaic(&scene(w, h), w, [0, 1, 1, 2], 0, 65535)); + spec.as_shot_neutral = Some([0.5, 1.0, 0.7]); + let mut s = Session::new(); + let r = open_bytes_as(&mut s, "IMG_0001.dng", &spec.build(), None, None).unwrap(); + assert!(r["warnings"].as_array().is_some_and(|w| w.iter().any(|m| m.as_str().is_some_and(|m| m.contains("DNG")))), "{r}"); + assert_eq!((doc(&s).size.width, doc(&s).size.height), (24, 16)); + assert_eq!(doc(&s).depth, photocraft_color::SampleType::U16); + // Damaged raw data is an error, not a crash. + let mut bad = spec.build(); + bad.truncate(bad.len() / 2); + assert!(open_bytes_as(&mut s, "bad.dng", &bad, None, None).is_err()); +} + #[test] fn revert_reloads_as_one_undoable_step() { let dir = tmp("revert"); diff --git a/crates/io/Cargo.toml b/crates/io/Cargo.toml index 6b9c1040..d18879ba 100644 --- a/crates/io/Cargo.toml +++ b/crates/io/Cargo.toml @@ -18,11 +18,13 @@ photocraft-doc = { workspace = true } photocraft-compose = { workspace = true } photocraft-psd = { workspace = true } photocraft-codecs = { workspace = true } +photocraft-raw = { workspace = true } photocraft-text = { workspace = true } photocraft-format = { path = "../format" } thiserror = { workspace = true } [dev-dependencies] photocraft-psd = { workspace = true, features = ["testgen"] } +photocraft-raw = { workspace = true, features = ["testgen"] } proptest = { workspace = true } photocraft-vector = { workspace = true } diff --git a/crates/io/src/flat.rs b/crates/io/src/flat.rs index 93ca736f..3372e5e0 100644 --- a/crates/io/src/flat.rs +++ b/crates/io/src/flat.rs @@ -13,6 +13,11 @@ use crate::{ExportOptions, ExportResult, ImportResult, IoError}; /// Decodes a flat image into a single-layer document. pub fn import_flat(name: &str, bytes: &[u8]) -> Result { let img = codecs::decode(bytes)?; + image_to_document(name, &img) +} + +/// A decoded flat image as a single-layer document. +pub(crate) fn image_to_document(name: &str, img: &Image) -> Result { let mut warnings = Vec::new(); let (mode, target_layout) = match img.layout() { ChannelLayout::Gray | ChannelLayout::GrayA => (ColorMode::Grayscale, ChannelLayout::GrayA), diff --git a/crates/io/src/lib.rs b/crates/io/src/lib.rs index 6f229f51..820dcc20 100644 --- a/crates/io/src/lib.rs +++ b/crates/io/src/lib.rs @@ -10,6 +10,9 @@ //! and smart-object data) and re-export writes them back; text, shape and //! smart-object layers also keep their pixels as the cached raster, and fill //! layers keep Photoshop's rendering in `Layer::fill_cache`. +//! * Camera raws (DNG, CR2, uncompressed / lossless TIFF-EP raws) via +//! `photocraft-raw`, developed into a 16-bit ProPhoto RGB "Background" +//! layer; unsupported raw variants fall back to the embedded JPEG preview. //! * Every other format goes through `photocraft-codecs` as a single //! "Background" layer (depth and Gray/RGB/CMYK model preserved). //! @@ -32,6 +35,7 @@ pub mod pattern_map; mod pixels; mod psd_export; mod psd_import; +pub mod raw; pub mod slices_map; pub mod text_styles_map; pub mod vector_map; @@ -63,6 +67,9 @@ pub enum IoError { /// Native `.pcraft` bundle failure. #[error("pcraft: {0}")] Pcraft(#[from] photocraft_format::FormatError), + /// Camera raw decode failure. + #[error("{0}")] + Raw(#[from] photocraft_raw::RawError), } /// Result of [`import`]. @@ -97,8 +104,8 @@ pub fn is_psd(bytes: &[u8]) -> bool { bytes.starts_with(b"8BPS") } -/// Imports a file. PSD/PSB are detected by magic; everything else is decoded -/// with `photocraft-codecs`. +/// Imports a file. PSD/PSB and camera raws are detected by magic; everything +/// else is decoded with `photocraft-codecs`. pub fn import(name: &str, bytes: &[u8]) -> Result { if photocraft_format::is_pcraft(bytes) { return Ok(ImportResult { document: photocraft_format::load_from_bytes(bytes)?, warnings: Vec::new() }); @@ -109,6 +116,9 @@ pub fn import(name: &str, bytes: &[u8]) -> Result { document.name = name.to_string(); return Ok(ImportResult { document, warnings }); } + if raw::is_raw(bytes) { + return raw::import_raw(name, bytes); + } flat::import_flat(name, bytes) } diff --git a/crates/io/src/raw.rs b/crates/io/src/raw.rs new file mode 100644 index 00000000..2220db77 --- /dev/null +++ b/crates/io/src/raw.rs @@ -0,0 +1,69 @@ +//! Camera raw files via `photocraft-raw`: developed to a 16-bit RGB document +//! in ProPhoto RGB (the embedded profile is the built-in ProPhoto-compatible +//! profile), or, for raw variants not decoded yet, the camera's embedded +//! JPEG preview with a warning. + +use std::sync::Arc; + +use photocraft_codecs::{self as codecs, ChannelLayout, Format, Image}; +use photocraft_raw::{DevelopOptions, Limits, RawError}; + +use crate::flat::image_to_document; +use crate::{ImportResult, IoError}; + +/// `true` if `bytes` are a camera raw file `photocraft-raw` recognises. +pub fn is_raw(bytes: &[u8]) -> bool { + photocraft_raw::is_raw(bytes) +} + +/// The decode limits shared with the flat codecs. +fn limits() -> Limits { + let l = codecs::Limits::default(); + Limits { max_width: l.max_width, max_height: l.max_height, max_pixels: l.max_pixels, max_alloc: l.max_alloc } +} + +/// Develops a raw file with the default settings. +pub fn import_raw(name: &str, bytes: &[u8]) -> Result { + import_raw_with(name, bytes, &DevelopOptions { limits: limits(), ..Default::default() }) +} + +/// Develops a raw file with explicit settings. +pub fn import_raw_with(name: &str, bytes: &[u8], opts: &DevelopOptions) -> Result { + let format = photocraft_raw::identify(bytes).map(|f| f.name()).unwrap_or("camera raw"); + match photocraft_raw::develop(bytes, opts) { + Ok(dev) => { + let img = Image::from_u16(dev.width, dev.height, ChannelLayout::Rgb, &dev.rgb)?; + let mut r = image_to_document(name, &img)?; + r.document.icc_profile = Some(Arc::new(photocraft_cms::Builtin::ProPhotoCompat.profile().to_bytes().to_vec())); + // "Canon" + "Canon EOS 80D" reads as "Canon EOS 80D". + let camera = match (dev.info.make.as_deref(), dev.info.model.as_deref()) { + (Some(make), Some(model)) if model.to_ascii_lowercase().starts_with(&make.to_ascii_lowercase()) => model.to_string(), + (make, model) => [make, model].into_iter().flatten().collect::>().join(" "), + }; + r.warnings.push(format!( + "{format}{} developed with default settings ({} demosaic, as-shot white balance) into 16-bit {}", + if camera.is_empty() { String::new() } else { format!(" from {camera}") }, + opts.demosaic.id(), + photocraft_raw::OUTPUT_SPACE + )); + r.warnings.extend(dev.warnings); + Ok(r) + } + Err(RawError::Unsupported(reason)) => match photocraft_raw::embedded_preview(bytes) { + Some(p) => { + let img = codecs::decode_as(Format::Jpeg, p.jpeg)?; + let mut r = image_to_document(name, &img)?; + r.warnings.insert( + 0, + format!( + "{format}: {reason} is not supported yet; opened the camera's embedded {}x{} JPEG preview instead (8-bit, not the raw sensor data)", + p.width, p.height + ), + ); + Ok(r) + } + None => Err(IoError::Raw(RawError::Unsupported(reason))), + }, + Err(e) => Err(IoError::Raw(e)), + } +} diff --git a/crates/io/tests/raw.rs b/crates/io/tests/raw.rs new file mode 100644 index 00000000..86cf1e77 --- /dev/null +++ b/crates/io/tests/raw.rs @@ -0,0 +1,90 @@ +//! Camera raw import: synthetic DNG / CR2 develop into a 16-bit ProPhoto +//! document; unsupported raw variants fall back to the embedded preview. + +use photocraft_codecs::{ChannelLayout, EncodeOptions, Format, Image}; +use photocraft_color::{ColorMode, SampleType}; +use photocraft_io::{IoError, import}; +use photocraft_raw::testgen::{Cr2Spec, DngSpec, TiffBuilder, Val, mosaic, scene}; + +#[test] +fn dng_opens_as_16_bit_prophoto() { + let (w, h) = (40, 24); + let mut spec = DngSpec::cfa(w, h, mosaic(&scene(w, h), w, [0, 1, 1, 2], 0, 65535)); + spec.as_shot_neutral = Some([0.6, 1.0, 0.8]); + let r = import("shot.dng", &spec.build()).unwrap(); + let d = &r.document; + assert_eq!((d.size.width, d.size.height), (40, 24)); + assert_eq!(d.mode, ColorMode::Rgb); + assert_eq!(d.depth, SampleType::U16); + assert_eq!(d.layers.len(), 1); + let icc = d.icc_profile.as_ref().expect("profile"); + assert_eq!(icc.as_slice(), &photocraft_cms::Builtin::ProPhotoCompat.profile().to_bytes()[..]); + assert!(r.warnings.iter().any(|w| w.contains("DNG") && w.contains("ProPhoto")), "{:?}", r.warnings); +} + +#[test] +fn cr2_opens() { + let (w, h) = (48, 16); + let data = mosaic(&scene(w, h), w, [0, 1, 1, 2], 256, 12000); + let spec = Cr2Spec { width: w, height: h, data, precision: 14, components: 2, slices: vec![24, 24], borders: None, wb_rggb: None, orientation: 6 }; + let r = import("IMG_0001.CR2", &spec.build()).unwrap(); + // Orientation 6 rotates the 48×16 sensor image to 16×48. + assert_eq!((r.document.size.width, r.document.size.height), (16, 48)); + assert_eq!(r.document.depth, SampleType::U16); +} + +/// A NEF-like file with Nikon's (undocumented) compression and a full-size +/// baseline JPEG preview in IFD0. +fn nef_with_preview() -> Vec { + let img = Image::from_u8(32, 20, ChannelLayout::Rgb, vec![180; 32 * 20 * 3]).unwrap(); + let jpeg = photocraft_codecs::encode(&img, Format::Jpeg, &EncodeOptions::default()).unwrap(); + let mut t = TiffBuilder::default(); + let strip = t.blob(vec![0; 64]); + let preview = t.blob(jpeg.clone()); + let raw = t.ifd(vec![ + (256, Val::Long(vec![8])), + (257, Val::Long(vec![8])), + (258, Val::Short(vec![12])), + (259, Val::Short(vec![34713])), + (262, Val::Short(vec![32803])), + (273, Val::Blobs(vec![strip])), + (279, Val::Long(vec![64])), + (33421, Val::Short(vec![2, 2])), + (33422, Val::Byte(vec![0, 1, 1, 2])), + ]); + let ifd0 = t.ifd(vec![ + (271, Val::Ascii("NIKON CORPORATION".into())), + (330, Val::Ifds(vec![raw])), + (513, Val::Blobs(vec![preview])), + (514, Val::Long(vec![jpeg.len() as u32])), + ]); + t.chain = vec![ifd0]; + t.build() +} + +#[test] +fn unsupported_raw_falls_back_to_the_embedded_preview() { + let r = import("DSC_0001.NEF", &nef_with_preview()).unwrap(); + assert_eq!((r.document.size.width, r.document.size.height), (32, 20)); + assert!(r.warnings.first().is_some_and(|w| w.contains("Nikon compressed NEF") && w.contains("embedded")), "{:?}", r.warnings); +} + +#[test] +fn unsupported_raw_without_preview_is_a_clear_error() { + let mut cr3 = vec![0, 0, 0, 24]; + cr3.extend_from_slice(b"ftypcrx "); + cr3.extend_from_slice(&[0; 12]); + match import("IMG_0001.CR3", &cr3) { + Err(e @ IoError::Raw(_)) => assert!(e.to_string().contains("CR3"), "{e}"), + Err(e) => panic!("expected a raw error, got {e}"), + Ok(_) => panic!("CR3 must not decode yet"), + } +} + +#[test] +fn truncated_raws_do_not_panic() { + let b = nef_with_preview(); + for n in (0..b.len()).step_by(7) { + let _ = import("x.nef", &b[..n]); + } +} diff --git a/crates/raw/Cargo.toml b/crates/raw/Cargo.toml new file mode 100644 index 00000000..90d4c83d --- /dev/null +++ b/crates/raw/Cargo.toml @@ -0,0 +1,27 @@ +[package] +name = "photocraft-raw" +description = "Clean-room, pure-Rust camera raw decoder and developer (DNG, CR2, TIFF/EP raws)" +version.workspace = true +edition.workspace = true +license.workspace = true +rust-version.workspace = true + +[features] +default = [] +# Synthetic raw-file generators (DNG / CR2 / TIFF-EP writers and a lossless-JPEG +# encoder) for tests and fuzz seeds. +testgen = [] + +[dependencies] + +# Parallel bands on native; wasm builds stay single-threaded. +[target.'cfg(not(target_arch = "wasm32"))'.dependencies] +rayon = { workspace = true } + +[dev-dependencies] +photocraft-raw = { path = ".", features = ["testgen"] } +# Writes sRGB PNG previews in the `rawinfo` example. +png = "0.18" + +[lints] +workspace = true diff --git a/crates/raw/README.md b/crates/raw/README.md new file mode 100644 index 00000000..f868a7d6 --- /dev/null +++ b/crates/raw/README.md @@ -0,0 +1,69 @@ +# photocraft-raw + +A clean-room, pure-Rust camera raw decoder and developer. The crate is standalone (no workspace +dependencies), has no `unsafe`, does no I/O (`&[u8]` in), builds for `wasm32-unknown-unknown` +(sequential there, rayon-parallel on native) and never panics on hostile input: every offset is +bounds-checked and sizes are checked against `Limits` before allocating. + +```rust +use photocraft_raw::{develop, DevelopOptions, Demosaic}; + +let dev = develop(&bytes, &DevelopOptions { demosaic: Demosaic::Ahd, ..Default::default() })?; +// dev.rgb: interleaved 16-bit RGB in ProPhoto RGB (ROMM primaries, D50, gamma 1.8) +// dev.warnings: anything approximated or not applied +``` + +`photocraft-io` uses it so opening a raw file yields a normal 16-bit RGB document tagged with the +built-in ProPhoto-compatible profile. + +## Sources (clean-room) + +Implemented only from public specifications, papers and observation of files: + +* TIFF 6.0, TIFF/EP (ISO 12234-2) and the Adobe DNG Specification 1.7. +* ITU-T T.81 (ISO 10918-1) Annex H: lossless JPEG, process 14 ("LJ92"). +* The published description of Canon's CR2 container (header, raw IFD, slice tag 0xC640). +* Publicly documented maker-note / private tags: Canon SensorInfo (0x00E0) and ColorData + (0x4001), Nikon WB_RBLevels (0x000C) and BlackLevel (0x003D), Sony BlackLevel (0x7310) and + WB_RGGBLevels (0x7313). +* Demosaicing: Malvar, He & Cutler (ICASSP 2004); Hirakawa & Parks, "Adaptive + homogeneity-directed demosaicing" (IEEE TIP 2005). +* McCamy's CCT approximation (1992); the Bradford chromatic adaptation transform. + +No code from dcraw, LibRaw, rawspeed, rawler, rawloader or darktable was read or used, and no +camera colour tables were copied. + +## Support matrix + +| Format | Status | +|---|---| +| DNG | Uncompressed (8–16 bit, packed or not) and lossless JPEG; strips and tiles; CFA (Bayer) and LinearRaw; LinearizationTable, BlackLevel (+ repeat, DeltaH/V), WhiteLevel, ActiveArea, DefaultCrop, ColorMatrix1/2, CameraCalibration, ForwardMatrix, AnalogBalance, AsShotNeutral / AsShotWhiteXY, BaselineExposure, Orientation, OpcodeList2 GainMap (lens shading) | +| DNG (lossy JPEG, JPEG XL, floating point; opcodes other than GainMap) | Unsupported / not applied (reported) | +| CR2 | Lossless JPEG with slices, borders and as-shot white balance from the maker note, black measured on the masked border | +| CR2 sRAW / mRAW | Unsupported | +| NEF / NRW, ARW, PEF and other TIFF/EP raws | Uncompressed and lossless-JPEG (incl. Sony lossless ARW) CFA data | +| Nikon compressed NEF, Sony compressed ARW, Pentax compressed PEF | Unsupported (vendor compression; `photocraft-io` opens the embedded JPEG preview instead) | +| CR3, RAF, ORF, RW2 | Recognised, unsupported (preview fallback where a preview is found) | +| X-Trans and other non-Bayer CFAs | Unsupported | + +## Development pipeline + +1. Linearization table, per-position black level, scale to the white level. +2. White balance (as shot; or grey-world when the file has none; or explicit multipliers), + normalized so the smallest multiplier is 1, then clip to 1 so blown highlights stay white. +3. Demosaic: `Bilinear`, `Mhc` (Malvar–He–Cutler) or `Ahd` (default). +4. Camera → XYZ (D50) per the DNG specification (ColorMatrix interpolated by the white's + correlated colour temperature, or ForwardMatrix), → linear ProPhoto; exposure + (BaselineExposure + user EV); gamma 1.8; 16 bits. +5. Orientation. + +Files without colour calibration (CR2, NEF, ARW…) use a documented neutral fallback: the +white-balanced camera channels are treated as linear sRGB primaries (colours are plausible but +less saturated than a calibrated profile; converting to DNG gives calibrated colour). No tone +curve is applied: the result is a scene-referred rendering, flatter than a camera JPEG. + +## Tools + +`cargo run --release -p photocraft-raw --example rawinfo -- [--dump] [--demosaic ahd] [--png DIR] FILE...` +prints what was decoded, times decode and develop, and can write sRGB PNG previews and the +embedded JPEG previews. diff --git a/crates/raw/examples/mutate.rs b/crates/raw/examples/mutate.rs new file mode 100644 index 00000000..51413995 --- /dev/null +++ b/crates/raw/examples/mutate.rs @@ -0,0 +1,58 @@ +//! Mutation stress test: corrupt real or synthetic raw files at random and +//! check that decoding and developing never panics. +//! +//! ```text +//! cargo run --release -p photocraft-raw --example mutate -- ITERATIONS FILE... +//! ``` + +use photocraft_raw::{Demosaic, DevelopOptions, Limits, decode, develop_sensor, embedded_preview, identify}; + +fn main() { + let mut args = std::env::args().skip(1); + let n: usize = args.next().and_then(|s| s.parse().ok()).unwrap_or(200); + let limits = Limits { max_width: 16384, max_height: 16384, max_pixels: 1 << 26, max_alloc: 1 << 30 }; + let mut seed = 0x5EED_u64; + let mut rng = move || { + seed = seed.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407); + seed >> 33 + }; + let mut panics = 0; + for f in args { + let Ok(orig) = std::fs::read(&f) else { continue }; + for it in 0..n { + let mut b = orig.clone(); + match rng() % 4 { + 0 => b.truncate(rng() as usize % b.len().max(1)), + 1 => { + // Header / IFD region. + let span = b.len().clamp(1, 1 << 16); + for _ in 0..1 + rng() % 8 { + let at = rng() as usize % span; + b[at] = rng() as u8; + } + } + _ => { + for _ in 0..1 + rng() % 16 { + let at = rng() as usize % b.len().max(1); + b[at] = rng() as u8; + } + } + } + let r = std::panic::catch_unwind(|| { + let _ = identify(&b); + let _ = embedded_preview(&b); + if let Ok(s) = decode(&b, &limits) { + let _ = develop_sensor(&s, &DevelopOptions { demosaic: Demosaic::Bilinear, limits, ..Default::default() }); + } + }); + if r.is_err() { + panics += 1; + let path = format!("{f}.panic-{it}"); + eprintln!("PANIC on mutation {it} of {f}; saved {path}"); + let _ = std::fs::write(path, &b); + } + } + println!("{f}: {n} mutations done"); + } + println!("panics: {panics}"); +} diff --git a/crates/raw/examples/rawinfo.rs b/crates/raw/examples/rawinfo.rs new file mode 100644 index 00000000..04e29780 --- /dev/null +++ b/crates/raw/examples/rawinfo.rs @@ -0,0 +1,239 @@ +//! Inspect and time raw decoding. +//! +//! ```text +//! cargo run --release -p photocraft-raw --example rawinfo -- [--dump] [--demosaic bilinear|mhc|ahd] [--ppm out.ppm] FILE... +//! ``` +//! +//! `--ppm` writes the developed image (ProPhoto, gamma 1.8, 16-bit) as a binary PPM; +//! `--png` writes an sRGB 8-bit preview at most 1600 pixels wide (for eyeballing). + +use std::time::Instant; + +use photocraft_raw::{Demosaic, DevelopOptions, decode, develop_sensor, dump_structure, embedded_preview, identify}; + +fn main() { + let mut args = std::env::args().skip(1); + let mut dump = false; + let mut ppm: Option = None; + let mut png: Option = None; + let mut opts = DevelopOptions::default(); + let mut files = Vec::new(); + while let Some(a) = args.next() { + match a.as_str() { + "--dump" => dump = true, + "--ppm" => ppm = args.next(), + "--png" => png = args.next(), + "--demosaic" => opts.demosaic = args.next().and_then(|s| Demosaic::from_id(&s)).unwrap_or_default(), + "--write-synthetic" => { + // `--write-synthetic DIR`: write 24 MP synthetic DNG (uncompressed, LJ92 tiles) and CR2 files. + use photocraft_raw::testgen::{Cr2Spec, DngSpec, DngStorage, mosaic, scene}; + let dir = args.next().unwrap_or_else(|| ".".into()); + let (w, h) = (6000, 4000); + let data = mosaic(&scene(w, h), w, [0, 1, 1, 2], 512, 15000); + let mut spec = DngSpec::cfa(w, h, data.clone()); + spec.bits = 14; + spec.black = vec![512]; + spec.white = 15000; + spec.as_shot_neutral = Some([0.5, 1.0, 0.7]); + spec.color_matrix1 = Some((21, [0.9, -0.3, -0.1, -0.4, 1.3, 0.1, -0.1, 0.2, 0.6])); + let _ = std::fs::write(format!("{dir}/synthetic-uncompressed.dng"), spec.build()); + spec.storage = DngStorage::Lj92Tiles { width: 256, height: 256 }; + let _ = std::fs::write(format!("{dir}/synthetic-lj92.dng"), spec.build()); + let cr2 = Cr2Spec { + width: w, + height: h, + data, + precision: 14, + components: 4, + slices: vec![w / 2, w - w / 2], + borders: None, + wb_rggb: None, + orientation: 1, + }; + let _ = std::fs::write(format!("{dir}/synthetic.cr2"), cr2.build()); + } + "--synthetic" => { + // `--synthetic 6000x4000`: time a generated 14-bit DNG (uncompressed and + // lossless-JPEG tiled) and a CR2 of that size. + let dims = args.next().unwrap_or_default(); + let (w, h) = dims.split_once('x').and_then(|(a, b)| Some((a.parse().ok()?, b.parse().ok()?))).unwrap_or((6000, 4000)); + synthetic(w, h, &opts); + } + _ => files.push(a), + } + } + for f in files { + let bytes = match std::fs::read(&f) { + Ok(b) => b, + Err(e) => { + eprintln!("{f}: {e}"); + continue; + } + }; + println!("== {f}: {:?}", identify(&bytes)); + if dump { + println!("{}", dump_structure(&bytes)); + } + if let Some(p) = embedded_preview(&bytes) { + println!(" preview: {}x{} JPEG, {} bytes", p.width, p.height, p.jpeg.len()); + if let Some(dir) = png.as_ref().filter(|d| std::path::Path::new(d).is_dir()) { + let stem = std::path::Path::new(&f).file_stem().map(|s| s.to_string_lossy().into_owned()).unwrap_or_default(); + let _ = std::fs::write(format!("{dir}/{stem}.preview.jpg"), p.jpeg); + } + } + let t0 = Instant::now(); + let sensor = match decode(&bytes, &opts.limits) { + Ok(s) => s, + Err(e) => { + println!(" decode: {e}"); + continue; + } + }; + let t1 = Instant::now(); + if std::env::var_os("RAWINFO_COLUMNS").is_some() { + // Mean of each of the first 200 columns over the middle half of the rows. + let (w, h) = (sensor.width, sensor.height); + let means: Vec = (0..200.min(w)) + .map(|x| { + let rows = h / 4..3 * h / 4; + let n = rows.len().max(1); + let s: u64 = rows.map(|y| u64::from(sensor.data[y * w + x])).sum(); + format!("{x}:{}", s / n as u64) + }) + .collect(); + println!(" columns {}", means.join(" ")); + } + println!( + " sensor {}x{}x{} cfa {:?} active {:?} crop {:?} black {:?} white {:?} wb {:?} neutral {:?} orient {} BE {}", + sensor.width, + sensor.height, + sensor.samples, + sensor.cfa.as_ref().map(|c| c.phase(sensor.crop.x, sensor.crop.y)), + sensor.active, + sensor.crop, + &sensor.black.values[..sensor.black.values.len().min(4)], + sensor.white, + sensor.camera_wb, + sensor.color.as_shot_neutral, + sensor.orientation, + sensor.baseline_exposure + ); + let dev = match develop_sensor(&sensor, &opts) { + Ok(d) => d, + Err(e) => { + println!(" develop: {e}"); + continue; + } + }; + let t2 = Instant::now(); + println!( + " developed {}x{} ({:?}) decode {:.0} ms, develop {:.0} ms; wb {:?}; warnings {:?}", + dev.width, + dev.height, + opts.demosaic, + (t1 - t0).as_secs_f64() * 1e3, + (t2 - t1).as_secs_f64() * 1e3, + dev.info.wb_multipliers, + dev.warnings + ); + if let Some(p) = &png { + // A directory gets one `.png` per input. + let path = if std::path::Path::new(p).is_dir() { + let stem = std::path::Path::new(&f).file_stem().map(|s| s.to_string_lossy().into_owned()).unwrap_or_default(); + format!("{p}/{stem}.png") + } else { + p.clone() + }; + write_png(&path, &dev); + } + if let Some(p) = &ppm { + let mut out = format!("P6\n{} {}\n65535\n", dev.width, dev.height).into_bytes(); + for v in &dev.rgb { + out.extend_from_slice(&v.to_be_bytes()); + } + if let Err(e) = std::fs::write(p, out) { + eprintln!("{p}: {e}"); + } + } + } +} + +fn synthetic(w: usize, h: usize, opts: &DevelopOptions) { + use photocraft_raw::testgen::{Cr2Spec, DngSpec, DngStorage, mosaic, scene}; + let data = mosaic(&scene(w, h), w, [0, 1, 1, 2], 512, 15000); + let mut spec = DngSpec::cfa(w, h, data.clone()); + spec.bits = 14; + spec.black = vec![512]; + spec.white = 15000; + spec.as_shot_neutral = Some([0.5, 1.0, 0.7]); + spec.color_matrix1 = Some((21, [0.9, -0.3, -0.1, -0.4, 1.3, 0.1, -0.1, 0.2, 0.6])); + let cr2 = Cr2Spec { width: w, height: h, data, precision: 14, components: 4, slices: vec![w / 2, w - w / 2], borders: None, wb_rggb: None, orientation: 1 }; + let files = [ + ("DNG uncompressed", spec.build()), + ("DNG LJ92 256x256 tiles", { + spec.storage = DngStorage::Lj92Tiles { width: 256, height: 256 }; + spec.build() + }), + ("CR2 (one LJ92 stream)", cr2.build()), + ]; + for (name, b) in files { + for m in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] { + let t0 = Instant::now(); + let s = decode(&b, &opts.limits).expect("synthetic decode"); + let t1 = Instant::now(); + let d = develop_sensor(&s, &DevelopOptions { demosaic: m, ..opts.clone() }).expect("synthetic develop"); + let t2 = Instant::now(); + println!( + "{name:24} {w}x{h} {:8}: decode {:6.0} ms + develop {:6.0} ms = {:6.0} ms ({} px out)", + m.id(), + (t1 - t0).as_secs_f64() * 1e3, + (t2 - t1).as_secs_f64() * 1e3, + (t2 - t0).as_secs_f64() * 1e3, + d.rgb.len() / 3 + ); + } + } +} + +/// ProPhoto (gamma 1.8) → sRGB 8-bit, box-downscaled to at most 1600 px wide. +fn write_png(path: &str, dev: &photocraft_raw::Developed) { + // Linear ProPhoto → XYZ D50 → linear sRGB (Bradford-adapted D50 → D65), combined. + const M: [[f32; 3]; 3] = [[2.0341, -0.7276, -0.3065], [-0.2289, 1.2317, -0.0028], [-0.0086, -0.1534, 1.1620]]; + let (w, h) = (dev.width as usize, dev.height as usize); + // RAWINFO_CROP="x,y,w,h" writes that region at 100 % instead of a downscaled whole. + let crop: Option> = std::env::var("RAWINFO_CROP").ok().map(|s| s.split(',').filter_map(|v| v.trim().parse().ok()).collect()); + let (x0, y0, k, ow, oh) = match crop.as_deref() { + Some(&[x, y, cw, ch]) if x + cw <= w && y + ch <= h => (x, y, 1, cw, ch), + _ => { + let k = w.div_ceil(1600).max(1); + (0, 0, k, w / k, h / k) + } + }; + let mut out = Vec::with_capacity(ow * oh * 3); + for oy in 0..oh { + for ox in 0..ow { + let mut acc = [0.0f32; 3]; + for dy in 0..k { + for dx in 0..k { + let i = ((y0 + oy * k + dy) * w + x0 + ox * k + dx) * 3; + for (c, a) in acc.iter_mut().enumerate() { + *a += (f32::from(dev.rgb[i + c]) / 65535.0).powf(1.8); + } + } + } + let lin = acc.map(|v| v / (k * k) as f32); + for row in M { + let v = (row[0] * lin[0] + row[1] * lin[1] + row[2] * lin[2]).clamp(0.0, 1.0); + let e = if v <= 0.0031308 { 12.92 * v } else { 1.055 * v.powf(1.0 / 2.4) - 0.055 }; + out.push((e * 255.0 + 0.5) as u8); + } + } + } + let Ok(f) = std::fs::File::create(path) else { return eprintln!("{path}: cannot create") }; + let mut enc = png::Encoder::new(std::io::BufWriter::new(f), ow as u32, oh as u32); + enc.set_color(png::ColorType::Rgb); + enc.set_depth(png::BitDepth::Eight); + if let Err(e) = enc.write_header().and_then(|mut w| w.write_image_data(&out)) { + eprintln!("{path}: {e}"); + } +} diff --git a/crates/raw/fuzz/Cargo.toml b/crates/raw/fuzz/Cargo.toml new file mode 100644 index 00000000..5917eb48 --- /dev/null +++ b/crates/raw/fuzz/Cargo.toml @@ -0,0 +1,34 @@ +[package] +name = "photocraft-raw-fuzz" +version = "0.0.0" +publish = false +edition = "2024" +license = "MIT OR Apache-2.0" + +[package.metadata] +cargo-fuzz = true + +# Standalone workspace so the main workspace's `crates/*` glob and lints do +# not apply (cargo-fuzz needs nightly + sanitizer flags). +[workspace] + +[dependencies] +libfuzzer-sys = "0.4" +photocraft-raw = { path = ".." } + +[profile.release] +debug = 1 + +[[bin]] +name = "develop" +path = "fuzz_targets/develop.rs" +test = false +doc = false +bench = false + +[[bin]] +name = "ljpeg" +path = "fuzz_targets/ljpeg.rs" +test = false +doc = false +bench = false diff --git a/crates/raw/fuzz/fuzz_targets/develop.rs b/crates/raw/fuzz/fuzz_targets/develop.rs new file mode 100644 index 00000000..31c2394e --- /dev/null +++ b/crates/raw/fuzz/fuzz_targets/develop.rs @@ -0,0 +1,15 @@ +#![no_main] + +use libfuzzer_sys::fuzz_target; +use photocraft_raw::{Demosaic, DevelopOptions, Limits, decode, develop_sensor, embedded_preview, identify}; + +// Any bytes: identify, find the preview, decode and develop (fast demosaic) +// under tight limits. Must never panic or exceed the limits. +fuzz_target!(|data: &[u8]| { + let limits = Limits { max_width: 4096, max_height: 4096, max_pixels: 1 << 20, max_alloc: 64 << 20 }; + let _ = identify(data); + let _ = embedded_preview(data); + if let Ok(s) = decode(data, &limits) { + let _ = develop_sensor(&s, &DevelopOptions { demosaic: Demosaic::Bilinear, limits, ..Default::default() }); + } +}); diff --git a/crates/raw/fuzz/fuzz_targets/ljpeg.rs b/crates/raw/fuzz/fuzz_targets/ljpeg.rs new file mode 100644 index 00000000..9f254ff3 --- /dev/null +++ b/crates/raw/fuzz/fuzz_targets/ljpeg.rs @@ -0,0 +1,8 @@ +#![no_main] + +use libfuzzer_sys::fuzz_target; + +// The lossless-JPEG decoder on its own: never panic, never allocate past the cap. +fuzz_target!(|data: &[u8]| { + let _ = photocraft_raw::decode_lossless_jpeg(data, 1 << 22); +}); diff --git a/crates/raw/src/color.rs b/crates/raw/src/color.rs new file mode 100644 index 00000000..dc6d977c --- /dev/null +++ b/crates/raw/src/color.rs @@ -0,0 +1,337 @@ +//! Colour science for raw development: 3×3 matrix helpers, chromaticity +//! conversions, correlated colour temperature, Bradford adaptation and the +//! camera → XYZ (D50) → output-space mapping of the DNG specification +//! (Adobe DNG Specification 1.7, chapter 6 "Mapping Camera Color Space to +//! CIE XYZ Space"). + +pub type Mat3 = [[f64; 3]; 3]; + +pub const IDENTITY: Mat3 = [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]; + +/// CIE D50 white (the ICC PCS illuminant), XYZ with Y = 1. +pub const D50: [f64; 3] = [0.9642, 1.0, 0.8249]; +/// D50 chromaticity. +pub const D50_XY: [f64; 2] = [0.3457, 0.3585]; +/// D65 chromaticity. +pub const D65_XY: [f64; 2] = [0.3127, 0.3290]; + +pub fn mul(a: &Mat3, b: &Mat3) -> Mat3 { + let mut o = [[0.0; 3]; 3]; + for (i, row) in o.iter_mut().enumerate() { + for (j, v) in row.iter_mut().enumerate() { + *v = (0..3).map(|k| a[i][k] * b[k][j]).sum(); + } + } + o +} + +pub fn apply(m: &Mat3, v: [f64; 3]) -> [f64; 3] { + [0, 1, 2].map(|i| m[i][0] * v[0] + m[i][1] * v[1] + m[i][2] * v[2]) +} + +pub fn diag(v: [f64; 3]) -> Mat3 { + [[v[0], 0.0, 0.0], [0.0, v[1], 0.0], [0.0, 0.0, v[2]]] +} + +pub fn scale(m: &Mat3, k: f64) -> Mat3 { + m.map(|r| r.map(|v| v * k)) +} + +/// Inverse, or `None` when (nearly) singular or not finite. +pub fn invert(m: &Mat3) -> Option { + let [[a, b, c], [d, e, f], [g, h, i]] = *m; + let det = a * (e * i - f * h) - b * (d * i - f * g) + c * (d * h - e * g); + if !det.is_finite() || det.abs() < 1e-12 { + return None; + } + let r = 1.0 / det; + let o = [ + [(e * i - f * h) * r, (c * h - b * i) * r, (b * f - c * e) * r], + [(f * g - d * i) * r, (a * i - c * g) * r, (c * d - a * f) * r], + [(d * h - e * g) * r, (b * g - a * h) * r, (a * e - b * d) * r], + ]; + o.iter().flatten().all(|v| v.is_finite()).then_some(o) +} + +/// A row-major 3×3 matrix from 9 values. +pub fn from_slice(v: &[f64]) -> Option { + if v.len() < 9 || v[..9].iter().any(|x| !x.is_finite()) { + return None; + } + Some([[v[0], v[1], v[2]], [v[3], v[4], v[5]], [v[6], v[7], v[8]]]) +} + +pub fn xy_to_xyz(xy: [f64; 2]) -> [f64; 3] { + let y = xy[1].max(1e-6); + [xy[0] / y, 1.0, (1.0 - xy[0] - xy[1]) / y] +} + +pub fn xyz_to_xy(v: [f64; 3]) -> Option<[f64; 2]> { + let s = v[0] + v[1] + v[2]; + (s.is_finite() && s > 1e-9).then(|| [v[0] / s, v[1] / s]) +} + +/// Correlated colour temperature from chromaticity (McCamy 1992, "Correlated +/// color temperature as an explicit function of chromaticity coordinates", +/// Color Res. Appl. 17(2)); clamped to 2000–50000 K. +pub fn cct(xy: [f64; 2]) -> f64 { + let n = (xy[0] - 0.3320) / (0.1858 - xy[1]); + let t = 449.0 * n.powi(3) + 3525.0 * n.powi(2) + 6823.3 * n + 5520.33; + if t.is_finite() { t.clamp(2000.0, 50000.0) } else { 5000.0 } +} + +/// Approximate temperature of an EXIF `LightSource` code (the +/// CalibrationIlluminant tags use the EXIF 2.3 table). `None` for unknown. +pub fn illuminant_temperature(code: u32) -> Option { + Some(match code { + 17 => 2856.0, // Standard light A + 18 => 4874.0, // Standard light B + 19 => 6774.0, // Standard light C + 20 => 5503.0, // D55 + 21 => 6504.0, // D65 + 22 => 7504.0, // D75 + 23 => 5003.0, // D50 + 24 => 3200.0, // ISO studio tungsten + 3 => 2850.0, // Tungsten + 1 | 4 | 9 => 5500.0, // Daylight, flash, fine weather + 10 => 6500.0, // Cloudy + 11 => 7500.0, // Shade + 12 => 6430.0, // Daylight fluorescent + 13 => 5000.0, // Day white fluorescent + 14 => 4150.0, // Cool white fluorescent + 15 => 3525.0, // White fluorescent + 16 => 2925.0, // Warm white fluorescent + 2 => 4200.0, // Fluorescent + _ => return None, + }) +} + +/// Bradford chromatic adaptation from white `src` to white `dst` (XYZ). +pub fn bradford(src: [f64; 3], dst: [f64; 3]) -> Mat3 { + const M: Mat3 = [[0.8951, 0.2664, -0.1614], [-0.7502, 1.7135, 0.0367], [0.0389, -0.0685, 1.0296]]; + let Some(mi) = invert(&M) else { return IDENTITY }; + let s = apply(&M, src); + let d = apply(&M, dst); + if s.iter().any(|v| v.abs() < 1e-9) { + return IDENTITY; + } + mul(&mi, &mul(&diag([d[0] / s[0], d[1] / s[1], d[2] / s[2]]), &M)) +} + +/// RGB → XYZ matrix for primaries and a white point (all as chromaticities), +/// white mapping to Y = 1. +pub fn rgb_to_xyz(prim: [[f64; 2]; 3], white: [f64; 2]) -> Mat3 { + let cols = prim.map(xy_to_xyz); + let p: Mat3 = [[cols[0][0], cols[1][0], cols[2][0]], [cols[0][1], cols[1][1], cols[2][1]], [cols[0][2], cols[1][2], cols[2][2]]]; + let Some(pi) = invert(&p) else { return IDENTITY }; + let s = apply(&pi, xy_to_xyz(white)); + mul(&p, &diag(s)) +} + +/// ROMM RGB (ProPhoto) primaries, ISO 22028-2; white D50. +pub const ROMM: [[f64; 2]; 3] = [[0.7347, 0.2653], [0.1596, 0.8404], [0.0366, 0.0001]]; +/// sRGB / Rec. 709 primaries. +pub const REC709: [[f64; 2]; 3] = [[0.64, 0.33], [0.30, 0.60], [0.15, 0.06]]; + +/// XYZ (D50) → linear ProPhoto RGB. +pub fn xyz_d50_to_prophoto() -> Mat3 { + invert(&rgb_to_xyz(ROMM, D50_XY)).unwrap_or(IDENTITY) +} + +/// Linear sRGB (D65) → XYZ adapted to D50: the documented neutral camera +/// matrix used when a file carries no colour calibration. +pub fn srgb_to_xyz_d50() -> Mat3 { + mul(&bradford(xy_to_xyz(D65_XY), xy_to_xyz(D50_XY)), &rgb_to_xyz(REC709, D65_XY)) +} + +/// One DNG calibration (ColorMatrixN with its companions). +#[derive(Debug, Clone, PartialEq)] +pub struct Calibration { + /// Illuminant temperature in kelvin. + pub temperature: f64, + /// XYZ → reference camera (ColorMatrixN). + pub color_matrix: Mat3, + /// White-balanced camera → XYZ D50 (ForwardMatrixN). + pub forward_matrix: Option, + /// Reference camera → individual camera (CameraCalibrationN). + pub camera_calibration: Mat3, +} + +/// Colour description of a raw file. +#[derive(Debug, Clone, PartialEq, Default)] +pub struct ColorInfo { + /// DNG calibrations (0, 1 or 2), sorted by temperature. + pub calibrations: Vec, + /// AnalogBalance. + pub analog_balance: Option<[f64; 3]>, + /// AsShotNeutral: camera coordinates of a neutral. + pub as_shot_neutral: Option<[f64; 3]>, + /// AsShotWhiteXY. + pub as_shot_white_xy: Option<[f64; 2]>, +} + +fn lerp_mat(a: &Mat3, b: &Mat3, g: f64) -> Mat3 { + let mut o = *a; + for i in 0..3 { + for j in 0..3 { + o[i][j] = a[i][j] * g + b[i][j] * (1.0 - g); + } + } + o +} + +impl ColorInfo { + /// Weight of the first calibration at temperature `t` (inverse-temperature + /// interpolation, DNG spec chapter 6). + fn weight(&self, t: f64) -> f64 { + match self.calibrations.as_slice() { + [a, b] if (a.temperature - b.temperature).abs() > 1.0 => { + let (t1, t2) = (a.temperature, b.temperature); + if t <= t1 { + 1.0 + } else if t >= t2 { + 0.0 + } else { + ((1.0 / t - 1.0 / t2) / (1.0 / t1 - 1.0 / t2)).clamp(0.0, 1.0) + } + } + _ => 1.0, + } + } + + /// The interpolated `(AB · CC · CM, ForwardMatrix, AB · CC)` at temperature `t`. + fn at(&self, t: f64) -> Option<(Mat3, Option, Mat3)> { + let first = self.calibrations.first()?; + let g = self.weight(t); + let pick = |f: &dyn Fn(&Calibration) -> Mat3| match self.calibrations.get(1) { + Some(b) => lerp_mat(&f(first), &f(b), g), + None => f(first), + }; + let cm = pick(&|c| c.color_matrix); + let cc = pick(&|c| c.camera_calibration); + let fm = match (first.forward_matrix, self.calibrations.get(1).map(|c| c.forward_matrix)) { + (Some(a), Some(Some(b))) => Some(lerp_mat(&a, &b, g)), + (Some(a), None) => Some(a), + _ => None, + }; + let ab = diag(self.analog_balance.unwrap_or([1.0; 3])); + let abcc = mul(&ab, &cc); + Some((mul(&abcc, &cm), fm, abcc)) + } + + /// Chromaticity of the white whose camera coordinates are `neutral`, + /// by fixed-point iteration over the temperature-dependent matrix. + pub fn neutral_to_xy(&self, neutral: [f64; 3]) -> [f64; 2] { + let mut xy = D50_XY; + for _ in 0..30 { + let Some((xyz_to_cam, _, _)) = self.at(cct(xy)) else { return xy }; + let Some(cam_to_xyz) = invert(&xyz_to_cam) else { return xy }; + let Some(next) = xyz_to_xy(apply(&cam_to_xyz, neutral)) else { return xy }; + if !(0.0..1.0).contains(&next[0]) || !(0.0..1.0).contains(&next[1]) { + return xy; + } + let done = (next[0] - xy[0]).abs() + (next[1] - xy[1]).abs() < 1e-7; + xy = next; + if done { + break; + } + } + xy + } + + /// Camera coordinates of a white with chromaticity `xy`, max component 1. + pub fn xy_to_neutral(&self, xy: [f64; 2]) -> Option<[f64; 3]> { + let (xyz_to_cam, _, _) = self.at(cct(xy))?; + let n = apply(&xyz_to_cam, xy_to_xyz(xy)); + let m = n.iter().cloned().fold(f64::MIN, f64::max); + (m > 1e-9 && n.iter().all(|v| v.is_finite() && *v > 0.0)).then(|| n.map(|v| v / m)) + } + + /// The matrix from white-balanced camera values (camera values divided by + /// `neutral / max(neutral)`) to XYZ D50, scaled so that the balanced + /// white (1, 1, 1) maps to D50 with Y = 1. `None` without calibrations. + pub fn balanced_to_xyz_d50(&self, neutral: [f64; 3]) -> Option { + let nmax = neutral.iter().cloned().fold(f64::MIN, f64::max); + if nmax.is_nan() || nmax <= 0.0 || neutral.iter().any(|v| v.is_nan() || *v <= 0.0) { + return None; + } + let xy = self.neutral_to_xy(neutral); + let (xyz_to_cam, fm, abcc) = self.at(cct(xy))?; + let unbalance = diag(neutral.map(|v| v / nmax)); + let m = match fm { + Some(fm) => { + // CameraToXYZ_D50 = FM · D · Inverse(AB · CC), D = Invert(diag(ReferenceNeutral)). + let abcc_inv = invert(&abcc)?; + let reference = apply(&abcc_inv, neutral); + if reference.iter().any(|v| v.is_nan() || *v <= 0.0) { + return None; + } + let d = diag(reference.map(|v| 1.0 / v)); + mul(&mul(&fm, &mul(&d, &abcc_inv)), &unbalance) + } + None => { + let cam_to_xyz = invert(&xyz_to_cam)?; + let adapt = bradford(xy_to_xyz(xy), D50); + mul(&mul(&adapt, &cam_to_xyz), &unbalance) + } + }; + // Normalize: balanced white → Y = 1. + let y = apply(&m, [1.0; 3])[1]; + (y.is_finite() && y > 1e-9).then(|| scale(&m, 1.0 / y)) + } +} + +#[cfg(test)] +mod tests { + use super::*; + + fn close(a: f64, b: f64, eps: f64) -> bool { + (a - b).abs() < eps + } + + #[test] + fn inverse_roundtrip() { + let m = srgb_to_xyz_d50(); + let i = invert(&m).unwrap(); + let p = mul(&m, &i); + for (r, row) in p.iter().enumerate() { + for (c, v) in row.iter().enumerate() { + assert!(close(*v, if r == c { 1.0 } else { 0.0 }, 1e-9)); + } + } + assert!(invert(&[[0.0; 3]; 3]).is_none()); + } + + #[test] + fn white_points_map_to_d50() { + let w = apply(&srgb_to_xyz_d50(), [1.0; 3]); + assert!(close(w[0], 0.9642, 2e-3) && close(w[1], 1.0, 1e-9) && close(w[2], 0.8249, 2e-3), "{w:?}"); + let p = apply(&xyz_d50_to_prophoto(), xy_to_xyz(D50_XY)); + assert!(p.iter().all(|v| close(*v, 1.0, 1e-9)), "{p:?}"); + } + + #[test] + fn mccamy_reference_points() { + assert!(close(cct(D65_XY), 6504.0, 15.0)); + assert!(close(cct([0.44757, 0.40745]), 2856.0, 15.0)); // illuminant A + } + + #[test] + fn dng_color_matrix_roundtrip() { + // A synthetic camera whose XYZ → camera matrix is the inverse of sRGB: with a + // D65 neutral, the balanced pipeline must reproduce sRGB → XYZ D50. + let cam_to_xyz = rgb_to_xyz(REC709, D65_XY); + let cm = invert(&cam_to_xyz).unwrap(); + let info = ColorInfo { + calibrations: vec![Calibration { temperature: 6504.0, color_matrix: cm, forward_matrix: None, camera_calibration: IDENTITY }], + ..Default::default() + }; + let neutral = info.xy_to_neutral(D65_XY).unwrap(); + assert!(neutral.iter().all(|v| close(*v, 1.0, 1e-3)), "{neutral:?}"); + let xy = info.neutral_to_xy(neutral); + assert!(close(xy[0], D65_XY[0], 1e-4) && close(xy[1], D65_XY[1], 1e-4)); + let m = info.balanced_to_xyz_d50(neutral).unwrap(); + let w = apply(&m, [1.0; 3]); + assert!(close(w[0], 0.9642, 3e-3) && close(w[1], 1.0, 1e-9) && close(w[2], 0.8249, 3e-3), "{w:?}"); + } +} diff --git a/crates/raw/src/cr2.rs b/crates/raw/src/cr2.rs new file mode 100644 index 00000000..bd1c7d28 --- /dev/null +++ b/crates/raw/src/cr2.rs @@ -0,0 +1,209 @@ +//! Canon CR2. +//! +//! CR2 is a TIFF file whose header carries `CR`, the major and minor version +//! and the offset of the raw IFD (IFD3) at bytes 8–15. The raw IFD stores one +//! lossless JPEG (T.81 process 14) whose decoded sample sequence fills the +//! sensor image in vertical slices: the private tag 0xC640 gives the slice +//! count, slice width and last slice width; slice 0 is filled row by row, +//! then slice 1, and so on. sRAW/mRAW (subsampled YCbCr) is not supported. +//! +//! Sensor geometry and the as-shot white balance come from the Canon maker +//! note (an IFD at the MakerNote offset, values relative to the TIFF header): +//! tag 0x00E0 (SensorInfo: sensor size and the image borders) and tag +//! 0x4001 (ColorData, whose as-shot RGGB levels sit at a version-dependent +//! position; candidates are validated before use). The black level is +//! measured on the masked border; the white level is the clipping point +//! found in the data. + +use crate::error::{RawError, Result}; +use crate::sensor::{BlackLevels, Cfa, Rect, Sensor}; +use crate::tiff::{Ifd, Tiff, tag}; +use crate::{Limits, RawFormat, ljpeg}; + +const CANON_SENSOR_INFO: u16 = 0x00E0; +const CANON_COLOR_DATA: u16 = 0x4001; + +/// The Canon maker note IFD, if present. +pub(crate) fn maker_note(t: &Tiff, ifds: &[Ifd]) -> Option { + let e = ifds.iter().find_map(|i| i.get(tag::MAKER_NOTE).copied())?; + t.ifd_at(e.at, 0) +} + +/// As-shot white-balance multipliers (R, G, B) from ColorData. +fn as_shot_wb(t: &Tiff, mn: &Ifd) -> Option<[f64; 3]> { + let e = mn.get(CANON_COLOR_DATA)?; + let v = t.uints(e); + // Word offsets of WB_RGGBLevelsAsShot used by the ColorData versions. + for off in [0x3F, 0x47, 0x19, 0x22] { + let Some(l) = v.get(off..off + 4) else { continue }; + let [r, g1, g2, b] = [l[0], l[1], l[2], l[3]].map(f64::from); + let g = (g1 + g2) / 2.0; + let plausible = (256.0..16384.0).contains(&g) && (g1 - g2).abs() <= g * 0.05 && (0.25..8.0).contains(&(r / g)) && (0.25..8.0).contains(&(b / g)); + if plausible { + return Some([r / g, 1.0, b / g]); + } + } + None +} + +/// The white (clipping) level: the largest value when a meaningful share of +/// samples piles up at it, else the full range of the declared precision. +pub(crate) fn clip_level(data: &[u16], precision: u32) -> f32 { + let full = ((1u32 << precision.min(16)) - 1) as f32; + let step = (data.len() / 2_000_000).max(1); + let mut max = 0u16; + for v in data.iter().step_by(step) { + max = max.max(*v); + } + if max == 0 { + return full; + } + let near = data.iter().step_by(step).filter(|&&v| u32::from(v) + u32::from(max / 256) >= u32::from(max)).count(); + let sampled = data.len().div_ceil(step); + if near * 2000 >= sampled { + // At least 0.05 % of samples at the top: that is the clipping point. + // Back off slightly so the clipped plateau maps to white. + f32::from(max) * 0.995 + } else { + full + } +} + +/// Mean of a region per 2×2 position. +pub(crate) fn masked_black(data: &[u16], width: usize, area: Rect) -> Option<[f32; 4]> { + if area.is_empty() || area.width < 2 || area.height < 2 { + return None; + } + let mut sum = [0u64; 4]; + let mut n = [0u64; 4]; + for y in area.y..area.y + area.height { + for x in area.x..area.x + area.width { + let v = *data.get(y * width + x)?; + let k = (y & 1) * 2 + (x & 1); + sum[k] += u64::from(v); + n[k] += 1; + } + } + if n.contains(&0) { + return None; + } + Some([0, 1, 2, 3].map(|k| sum[k] as f32 / n[k] as f32)) +} + +pub(crate) fn decode(t: &Tiff, limits: &Limits) -> Result { + let ifds = t.all_ifds(); + let ifd0 = t.ifd_at(t.first_ifd, 0).ok_or_else(|| RawError::malformed("CR2 has no IFD0"))?; + let raw_off = t.u32_at(12).unwrap_or(0) as usize; + let raw = t.ifd_at(raw_off, 0).ok_or_else(|| RawError::malformed("CR2 raw IFD is missing"))?; + let off = t.tag_uint(&raw, tag::STRIP_OFFSETS).ok_or_else(|| RawError::malformed("CR2 raw data offset missing"))? as usize; + let len = t.tag_uint(&raw, tag::STRIP_BYTE_COUNTS).map(|l| l as usize).unwrap_or(t.data.len().saturating_sub(off)); + let src = t.bytes(off, len).or_else(|| t.data.get(off..)).ok_or_else(|| RawError::malformed("CR2 raw data lies outside the file"))?; + + let frame = ljpeg::frame(src).map_err(|e| match e { + RawError::Unsupported(m) if m.contains("subsampled") => RawError::unsupported("Canon sRAW / mRAW"), + e => e, + })?; + let total = frame.samples().ok_or_else(|| RawError::malformed("CR2 frame too large"))?; + let slices = t.tag_uints(&raw, tag::CR2_SLICE); + let (width, widths) = match slices.as_slice() { + [n, w, last] if *n < 64 && *w > 0 && *last > 0 => { + let (n, w, last) = (*n as usize, *w as usize, *last as usize); + let mut ws = vec![w; n]; + ws.push(last); + (n * w + last, ws) + } + _ => { + let w = frame.width * frame.components; + (w, vec![w]) + } + }; + if width == 0 || total % width != 0 { + return Err(RawError::malformed("CR2 slices do not match the JPEG frame")); + } + let height = total / width; + limits.check(width as u64, height as u64, 2)?; + let (_, samples) = ljpeg::decode(src, total)?; + + // De-slice. + let mut data = vec![0u16; width * height]; + let mut pos = 0usize; + let mut x0 = 0usize; + for sw in widths { + for y in 0..height { + let src = samples.get(pos..pos + sw).ok_or_else(|| RawError::malformed("CR2 slice data is short"))?; + data[y * width + x0..y * width + x0 + sw].copy_from_slice(src); + pos += sw; + } + x0 += sw; + } + + let mut warnings = Vec::new(); + let mn = maker_note(t, &ifds); + let full = Rect::new(0, 0, width, height); + // SensorInfo: [1] width, [2] height, [5] left, [6] top, [7] right, [8] bottom border (inclusive). + let info = mn.as_ref().map(|m| t.tag_uints(m, CANON_SENSOR_INFO)).unwrap_or_default(); + let active = match info.get(5..9) { + Some(&[l, tp, r, b]) if r > l && b > tp => Rect::new(l as usize, tp as usize, (r - l + 1) as usize, (b - tp + 1) as usize).intersect(&full), + _ => full, + }; + let active = if active.is_empty() { full } else { active }; + // Black: the masked columns left of the image (skipping a few edge columns). + let black = if active.x >= 16 { + // Only the outer half: the columns next to the image area can be partly exposed. + let area = Rect::new(4, active.y, active.x / 2 - 4, active.height); + masked_black(&data, width, area) + } else { + None + }; + let black = match black { + Some(b) => BlackLevels { rows: 2, cols: 2, values: b.to_vec(), delta_h: Vec::new(), delta_v: Vec::new() }, + None => { + warnings.push("no masked sensor area; black level assumed to be 0".to_string()); + BlackLevels::uniform(0.0) + } + }; + let white = clip_level(&data, u32::from(frame.precision)); + // Canon sensors are RGGB from the top-left of the image area. + let cfa = Cfa { width: 2, height: 2, colors: vec![0, 1, 1, 2], origin_x: active.x, origin_y: active.y }; + let camera_wb = mn.as_ref().and_then(|m| as_shot_wb(t, m)); + // The black-level values index from the active area origin, but the masked + // measurement used data parity: realign them to active-area parity. + let black = realign_black(black, active); + Ok(Sensor { + format: RawFormat::Cr2, + make: t.tag_ascii(&ifd0, tag::MAKE), + model: t.tag_ascii(&ifd0, tag::MODEL), + width, + height, + samples: 1, + data, + cfa: Some(cfa), + linearization: None, + black, + white: [white; 3], + active, + crop: active, + color: Default::default(), + camera_wb, + orientation: t.tag_uint(&ifd0, tag::ORIENTATION).map(|o| o as u16).filter(|o| (1..=8).contains(o)).unwrap_or(1), + baseline_exposure: 0.0, + gain_maps: Vec::new(), + warnings, + }) +} + +/// Black levels measured per data-parity 2×2 position, re-indexed from the active origin. +pub(crate) fn realign_black(b: BlackLevels, active: Rect) -> BlackLevels { + if b.rows != 2 || b.cols != 2 || b.values.len() != 4 { + return b; + } + let v = &b.values; + let values = (0..4) + .map(|k| { + let (ax, ay) = (k & 1, k >> 1); + let (x, y) = ((active.x + ax) & 1, (active.y + ay) & 1); + v[y * 2 + x] + }) + .collect(); + BlackLevels { values, ..b } +} diff --git a/crates/raw/src/demosaic.rs b/crates/raw/src/demosaic.rs new file mode 100644 index 00000000..c6b7d332 --- /dev/null +++ b/crates/raw/src/demosaic.rs @@ -0,0 +1,500 @@ +//! Bayer demosaicing. +//! +//! * **Bilinear**: each missing colour is the mean of the same-colour +//! neighbours in the 3×3 window. +//! * **Malvar–He–Cutler**: H. S. Malvar, L.-W. He, R. Cutler, "High-quality +//! linear interpolation for demosaicing of Bayer-pattern color images", +//! ICASSP 2004: bilinear plus a gradient-correction term (5×5 kernels). +//! * **AHD**: K. Hirakawa, T. W. Parks, "Adaptive homogeneity-directed +//! demosaicing algorithm", IEEE Trans. Image Processing 14(3), 2005: +//! horizontally and vertically interpolated candidates, chosen per pixel by +//! their local homogeneity in CIELAB. +//! +//! Input is a [`Padded`] plane of normalized CFA samples (mirror-reflected +//! borders, which keep the CFA parity, so every kernel reads in bounds); +//! `phase` gives the colours (0 = R, 1 = G, 2 = B) of the 2×2 cell at the +//! image origin, row major. Output is interleaved RGB. Work is split in row +//! bands run in parallel. + +use crate::color::Mat3; +use crate::par; + +/// Demosaicing algorithm. +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum Demosaic { + /// Fastest; soft, with colour fringes on fine detail. + Bilinear, + /// Malvar–He–Cutler gradient-corrected linear interpolation: fast and sharp. + Mhc, + /// Adaptive homogeneity-directed (Hirakawa–Parks): best edges, slower. + #[default] + Ahd, +} + +impl Demosaic { + pub fn id(self) -> &'static str { + match self { + Demosaic::Bilinear => "bilinear", + Demosaic::Mhc => "mhc", + Demosaic::Ahd => "ahd", + } + } + + pub fn from_id(s: &str) -> Option { + match s.to_ascii_lowercase().as_str() { + "bilinear" | "fast" => Some(Demosaic::Bilinear), + "mhc" | "malvar" => Some(Demosaic::Mhc), + "ahd" | "best" => Some(Demosaic::Ahd), + _ => None, + } + } +} + +/// Border width of a [`Padded`] plane (even, so the CFA parity is unchanged). +pub(crate) const PAD: usize = 6; + +/// Mirror-reflects `i` into `0..n` (n ≥ 1); clamps when `n` is tiny. +#[inline] +fn reflect(i: isize, n: usize) -> usize { + let n = n as isize; + if n <= 1 { + return 0; + } + let mut i = i; + if i < 0 { + i = -i; + } + if i >= n { + i = 2 * (n - 1) - i; + } + i.clamp(0, n - 1) as usize +} + +/// A `w × h` plane stored with [`PAD`] extra pixels on every side. +pub(crate) struct Padded { + pub data: Vec, + pub w: usize, + pub h: usize, + pub stride: usize, +} + +impl Padded { + pub fn new(w: usize, h: usize) -> Self { + let stride = w + 2 * PAD; + Padded { data: vec![0.0; stride * (h + 2 * PAD)], w, h, stride } + } + + /// From an unpadded plane (tests and small inputs). + #[cfg(test)] + pub fn from_plane(plane: &[f32], w: usize, h: usize) -> Self { + let mut p = Padded::new(w, h); + for y in 0..h { + let at = (y + PAD) * p.stride + PAD; + p.data[at..at + w].copy_from_slice(&plane[y * w..(y + 1) * w]); + } + p.fill_borders(); + p + } + + /// Mutable interior rows, `rows` at a time, for parallel filling: calls + /// `f(first_row, rows_slice)` where each row is `stride` long and the + /// image pixels start at offset [`PAD`]. + pub fn fill_rows(&mut self, rows: usize, f: impl Fn(usize, &mut [f32]) + Sync + Send) { + let (s, h) = (self.stride, self.h); + let interior = &mut self.data[PAD * s..(PAD + h) * s]; + par::chunks_mut(interior, rows.max(1) * s, |b, chunk| f(b * rows.max(1), chunk)); + } + + /// Fills the borders by mirror reflection of the interior. + pub fn fill_borders(&mut self) { + let (w, h, s) = (self.w, self.h, self.stride); + if w == 0 || h == 0 { + return; + } + for y in PAD..PAD + h { + let row = &mut self.data[y * s..(y + 1) * s]; + for px in (0..PAD).chain(PAD + w..s) { + row[px] = row[PAD + reflect(px as isize - PAD as isize, w)]; + } + } + for py in (0..PAD).chain(PAD + h..h + 2 * PAD) { + let src = PAD + reflect(py as isize - PAD as isize, h); + self.data.copy_within(src * s..(src + 1) * s, py * s); + } + } +} + +/// Per 2×2 position, per output colour: kernel taps as (offset, weight). +fn kernels(phase: [u8; 4], method: Demosaic, stride: usize) -> [[Vec<(isize, f32)>; 3]; 4] { + let color = |x: isize, y: isize| phase[((y & 1) * 2 + (x & 1)) as usize]; + let s = stride as isize; + std::array::from_fn(|pos| { + let (px, py) = ((pos & 1) as isize, (pos >> 1) as isize); + let own = color(px, py); + std::array::from_fn(|c| { + let c = c as u8; + let taps: Vec<(isize, isize, f32)> = if c == own { + vec![(0, 0, 1.0)] + } else if method == Demosaic::Bilinear { + let t: Vec<(isize, isize)> = + (-1..=1).flat_map(|dy| (-1..=1).map(move |dx| (dx, dy))).filter(|&(dx, dy)| color(px + dx, py + dy) == c).collect(); + let w = 1.0 / t.len().max(1) as f32; + t.into_iter().map(|(dx, dy)| (dx, dy, w)).collect() + } else if c == 1 { + // Malvar–He–Cutler (weights / 8): G at R or B. + vec![(0, 0, 4.0), (-1, 0, 2.0), (1, 0, 2.0), (0, -1, 2.0), (0, 1, 2.0), (-2, 0, -1.0), (2, 0, -1.0), (0, -2, -1.0), (0, 2, -1.0)] + } else if own == 1 { + // R or B at G: from the row neighbours or the column neighbours. + let (ax, ay) = if color(px + 1, py) == c { (1, 0) } else { (0, 1) }; + vec![ + (0, 0, 5.0), + (-ax, -ay, 4.0), + (ax, ay, 4.0), + (-1, -1, -1.0), + (1, -1, -1.0), + (-1, 1, -1.0), + (1, 1, -1.0), + (-2 * ax, -2 * ay, -1.0), + (2 * ax, 2 * ay, -1.0), + (-2 * ay, -2 * ax, 0.5), + (2 * ay, 2 * ax, 0.5), + ] + } else { + // R at B or B at R. + vec![(0, 0, 6.0), (-1, -1, 2.0), (1, -1, 2.0), (-1, 1, 2.0), (1, 1, 2.0), (-2, 0, -1.5), (2, 0, -1.5), (0, -2, -1.5), (0, 2, -1.5)] + }; + let norm = if method == Demosaic::Bilinear || c == own { 1.0 } else { 1.0 / 8.0 }; + taps.into_iter().map(|(dx, dy, w)| (dy * s + dx, w * norm)).collect() + }) + }) +} + +/// Demosaics a padded CFA plane into interleaved RGB (`w × h × 3`). +pub(crate) fn demosaic(p: &Padded, phase: [u8; 4], method: Demosaic, to_xyz: &Mat3) -> Vec { + let (w, h, s) = (p.w, p.h, p.stride); + let mut out = vec![0.0f32; w * h * 3]; + if w == 0 || h == 0 || p.data.len() < s * (h + 2 * PAD) { + return out; + } + if method == Demosaic::Ahd { + ahd(p, phase, to_xyz, &mut out); + return out; + } + let ks = kernels(phase, method, s); + let rows = par::band_rows(w); + let d = &p.data; + par::chunks_mut(&mut out, rows * w * 3, |band, chunk| { + let y0 = band * rows; + for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() { + let y = y0 + r; + let base = (y + PAD) * s + PAD; + for (x, o) in row.as_chunks_mut::<3>().0.iter_mut().enumerate() { + let pos = (y & 1) * 2 + (x & 1); + let i = base + x; + if method == Demosaic::Bilinear { + // Direct Bayer bilinear (the kernels' result without the tap loops). + let own = usize::from(phase[pos]).min(2); + let (n, so, e, wv) = (d[i - s], d[i + s], d[i + 1], d[i - 1]); + if own == 1 { + let hc = usize::from(phase[(y & 1) * 2 + ((x + 1) & 1)]).min(2); + o[1] = d[i]; + o[hc] = (e + wv) * 0.5; + o[2 - hc] = (n + so) * 0.5; + } else { + o[own] = d[i]; + o[1] = (n + so + e + wv) * 0.25; + o[2 - own] = (d[i - s - 1] + d[i - s + 1] + d[i + s - 1] + d[i + s + 1]) * 0.25; + } + continue; + } + for (c, k) in ks[pos].iter().enumerate() { + let mut v = 0.0f32; + for &(off, wt) in k { + v += wt * d[(i as isize + off) as usize]; + } + o[c] = v.clamp(0.0, 1.0); + } + } + } + }); + out +} + +/// CIE f(t) for L*a*b*, tabulated on 0..2 with linear interpolation. +struct LabF { + table: Vec, +} + +const LABF_STEPS: usize = 1 << 14; +const LABF_RANGE: f32 = 2.0; + +impl LabF { + fn new() -> Self { + LabF { table: (0..=LABF_STEPS + 1).map(|i| Self::exact(i as f32 * LABF_RANGE / LABF_STEPS as f32)).collect() } + } + + fn exact(t: f32) -> f32 { + if t > 0.008856 { t.cbrt() } else { 7.787 * t + 16.0 / 116.0 } + } + + #[inline] + fn f(&self, t: f32) -> f32 { + if !(0.0..LABF_RANGE).contains(&t) { + return Self::exact(t.max(0.0)); + } + let x = t * (LABF_STEPS as f32 / LABF_RANGE); + let i = (x as usize).min(LABF_STEPS); + let a = self.table[i]; + a + (self.table[i + 1] - a) * (x - i as f32) + } +} + +/// Output rows per AHD band. +const AHD_BAND: usize = 32; +/// Extra band rows above and below the output rows. +const M: usize = 3; + +fn ahd(p: &Padded, phase: [u8; 4], to_xyz: &Mat3, out: &mut [f32]) { + let (w, s) = (p.w, p.stride); + let d = &p.data; + // Colour at padded coordinates (PAD is even, so the phase is unchanged). + let color = |px: usize, py: usize| phase[(py & 1) * 2 + (px & 1)]; + let white = [0.9642f32, 1.0, 0.8249]; + let m: [[f32; 3]; 3] = std::array::from_fn(|i| std::array::from_fn(|j| to_xyz[i][j] as f32 / white[i])); + let labf = LabF::new(); + par::chunks_mut(out, AHD_BAND * w * 3, |band, chunk| { + let rows = chunk.len() / (w * 3); + let n = rows + 2 * M; + // Band row b ↔ padded row py(b); padded columns are used as is. + let py = |b: usize| band * AHD_BAND + PAD - M + b; + // 1. Green interpolated horizontally (0) and vertically (1), columns 2..s-2. + let mut g = [vec![0.0f32; n * s], vec![0.0f32; n * s]]; + for b in 0..n { + let y = py(b); + for x in 2..s - 2 { + let i = y * s + x; + let c = d[i]; + let k = b * s + x; + if color(x, y) == 1 { + g[0][k] = c; + g[1][k] = c; + continue; + } + let (l, r) = (d[i - 1], d[i + 1]); + let gh = (l + r) * 0.5 + (2.0 * c - d[i - 2] - d[i + 2]) * 0.25; + g[0][k] = gh.clamp(l.min(r), l.max(r)); + let (u, dn) = (d[i - s], d[i + s]); + let gv = (u + dn) * 0.5 + (2.0 * c - d[i - 2 * s] - d[i + 2 * s]) * 0.25; + g[1][k] = gv.clamp(u.min(dn), u.max(dn)); + } + } + // 2. RGB candidates and their Lab, band rows 1..n-1, columns 3..s-3. + let mut rgb = [vec![0.0f32; n * s * 3], vec![0.0f32; n * s * 3]]; + let mut lab = [vec![0.0f32; n * s * 3], vec![0.0f32; n * s * 3]]; + for dir in 0..2 { + let gd = &g[dir]; + // Colour difference (sample − green) at band row b, column x. + let diff = |b: usize, x: usize| d[py(b) * s + x] - gd[b * s + x]; + let (rgb_d, lab_d) = (&mut rgb[dir], &mut lab[dir]); + for b in 1..n - 1 { + let y = py(b); + for x in 3..s - 3 { + let own = color(x, y); + let gv = gd[b * s + x]; + let mut px = [gv; 3]; + if own == 1 { + let hc = usize::from(color(x + 1, y)).min(2); + px[hc] = gv + 0.5 * (diff(b, x - 1) + diff(b, x + 1)); + px[2 - hc] = gv + 0.5 * (diff(b - 1, x) + diff(b + 1, x)); + } else { + let o = usize::from(own).min(2); + px[o] = d[y * s + x]; + px[2 - o] = gv + 0.25 * (diff(b - 1, x - 1) + diff(b - 1, x + 1) + diff(b + 1, x - 1) + diff(b + 1, x + 1)); + } + let px = px.map(|v| v.clamp(0.0, 1.0)); + let i = (b * s + x) * 3; + rgb_d[i..i + 3].copy_from_slice(&px); + let xyz: [f32; 3] = std::array::from_fn(|k| m[k][0] * px[0] + m[k][1] * px[1] + m[k][2] * px[2]); + let (fx, fy, fz) = (labf.f(xyz[0]), labf.f(xyz[1]), labf.f(xyz[2])); + lab_d[i] = 116.0 * fy - 16.0; + lab_d[i + 1] = 500.0 * (fx - fy); + lab_d[i + 2] = 200.0 * (fy - fz); + } + } + } + // 3. Homogeneity maps, band rows 2..n-2, columns 4..s-4. + let mut hom = [vec![0u8; n * s], vec![0u8; n * s]]; + let at = |l: &[f32], k: usize| (l[k * 3], l[k * 3 + 1], l[k * 3 + 2]); + let dist = |a: (f32, f32, f32), c: (f32, f32, f32)| ((a.0 - c.0).abs(), (a.1 - c.1) * (a.1 - c.1) + (a.2 - c.2) * (a.2 - c.2)); + for b in 2..n - 2 { + for x in 4..s - 4 { + let k = b * s + x; + let nb = [k - 1, k + 1, k - s, k + s]; + let (ch, cv) = (at(&lab[0], k), at(&lab[1], k)); + let dh: [(f32, f32); 4] = nb.map(|j| dist(ch, at(&lab[0], j))); + let dv: [(f32, f32); 4] = nb.map(|j| dist(cv, at(&lab[1], j))); + // The horizontal candidate judged along the row, the vertical along the column. + let eps_l = dh[0].0.max(dh[1].0).min(dv[2].0.max(dv[3].0)); + let eps_c = dh[0].1.max(dh[1].1).min(dv[2].1.max(dv[3].1)); + hom[0][k] = dh.iter().filter(|(l, c)| *l <= eps_l && *c <= eps_c).count() as u8; + hom[1][k] = dv.iter().filter(|(l, c)| *l <= eps_l && *c <= eps_c).count() as u8; + } + } + // 4. Per pixel, the candidate with the larger 3×3 homogeneity sum. + for r in 0..rows { + let b = r + M; + for x in 0..w { + let px = x + PAD; + let (mut sh, mut sv) = (0u32, 0u32); + for bb in b - 1..=b + 1 { + for xx in px - 1..=px + 1 { + sh += u32::from(hom[0][bb * s + xx]); + sv += u32::from(hom[1][bb * s + xx]); + } + } + let i = (b * s + px) * 3; + let o = &mut chunk[(r * w + x) * 3..(r * w + x) * 3 + 3]; + for k in 0..3 { + o[k] = match sh.cmp(&sv) { + std::cmp::Ordering::Greater => rgb[0][i + k], + std::cmp::Ordering::Less => rgb[1][i + k], + std::cmp::Ordering::Equal => 0.5 * (rgb[0][i + k] + rgb[1][i + k]), + }; + } + } + } + }); +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::color::IDENTITY; + + /// Mosaics an RGB image with the given phase. + fn mosaic(rgb: &[[f32; 3]], w: usize, h: usize, phase: [u8; 4]) -> Vec { + (0..w * h).map(|i| rgb[i][phase[((i / w) & 1) * 2 + ((i % w) & 1)] as usize]).collect() + } + + fn run(cfa: &[f32], w: usize, h: usize, phase: [u8; 4], m: Demosaic) -> Vec { + demosaic(&Padded::from_plane(cfa, w, h), phase, m, &IDENTITY) + } + + fn psnr(a: &[f32], b: &[[f32; 3]], w: usize, h: usize, border: usize) -> f64 { + let mut se = 0.0f64; + let mut n = 0usize; + for y in border..h - border { + for x in border..w - border { + for c in 0..3 { + let d = f64::from(a[(y * w + x) * 3 + c] - b[y * w + x][c]); + se += d * d; + n += 1; + } + } + } + 10.0 * (1.0 / (se / n as f64).max(1e-12)).log10() + } + + #[test] + fn padding_reflects() { + let p = Padded::from_plane(&[1.0, 2.0, 3.0, 4.0, 5.0, 6.0], 3, 2); + let row = |y: usize| p.data[(y + PAD) * p.stride + PAD - 2..(y + PAD) * p.stride + PAD + 5].to_vec(); + assert_eq!(row(0), vec![3.0, 2.0, 1.0, 2.0, 3.0, 2.0, 1.0]); + // Row −1 mirrors row 1. + let r = PAD - 1; + assert_eq!(p.data[r * p.stride + PAD..r * p.stride + PAD + 3].to_vec(), vec![4.0, 5.0, 6.0]); + } + + #[test] + fn flat_colour_is_exact() { + let (w, h) = (16, 12); + let rgb = vec![[0.2f32, 0.5, 0.8]; w * h]; + for phase in [[0, 1, 1, 2], [2, 1, 1, 0], [1, 0, 2, 1], [1, 2, 0, 1]] { + for m in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] { + let out = run(&mosaic(&rgb, w, h, phase), w, h, phase, m); + for (i, v) in out.as_chunks::<3>().0.iter().enumerate() { + for c in 0..3 { + assert!((v[c] - rgb[i][c]).abs() < 1e-5, "{m:?} {phase:?} px {i} c {c}: {}", v[c]); + } + } + } + } + } + + #[test] + fn smooth_gradient_quality() { + // A smooth colour field: every method should be accurate. + let (w, h) = (64, 48); + let rgb: Vec<[f32; 3]> = (0..w * h) + .map(|i| { + let (x, y) = ((i % w) as f32 / w as f32, (i / w) as f32 / h as f32); + [0.2 + 0.6 * x, 0.3 + 0.4 * (x * y), 0.8 - 0.5 * y] + }) + .collect(); + let phase = [0, 1, 1, 2]; + let cfa = mosaic(&rgb, w, h, phase); + let p: Vec = [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd].iter().map(|&m| psnr(&run(&cfa, w, h, phase, m), &rgb, w, h, 2)).collect(); + assert!(p.iter().all(|&v| v > 40.0), "{p:?}"); + } + + #[test] + fn edges_ahd_beats_bilinear() { + // Vertical and horizontal black/white stripes of width 3: AHD's + // directional choice must reduce error versus bilinear. + let (w, h) = (48, 48); + let rgb: Vec<[f32; 3]> = (0..w * h) + .map(|i| { + let (x, y) = (i % w, i / w); + let v = if y < h / 2 { ((x / 3) % 2) as f32 } else { ((y / 3) % 2) as f32 }; + [0.1 + 0.8 * v, 0.1 + 0.8 * v, 0.1 + 0.8 * v] + }) + .collect(); + let phase = [0, 1, 1, 2]; + let cfa = mosaic(&rgb, w, h, phase); + let bil = psnr(&run(&cfa, w, h, phase, Demosaic::Bilinear), &rgb, w, h, 3); + let ahd = psnr(&run(&cfa, w, h, phase, Demosaic::Ahd), &rgb, w, h, 3); + let mhc = psnr(&run(&cfa, w, h, phase, Demosaic::Mhc), &rgb, w, h, 3); + assert!(ahd > bil + 3.0, "ahd {ahd:.1} dB vs bilinear {bil:.1} dB"); + assert!(mhc > bil, "mhc {mhc:.1} dB vs bilinear {bil:.1} dB"); + } + + #[test] + fn bilinear_fast_path_matches_kernels() { + // The direct bilinear code must equal the generic kernel evaluation. + let (w, h) = (20, 14); + let cfa: Vec = (0..w * h).map(|i| ((i * 7919) % 1000) as f32 / 1000.0).collect(); + for phase in [[0, 1, 1, 2], [1, 2, 0, 1]] { + let p = Padded::from_plane(&cfa, w, h); + let fast = demosaic(&p, phase, Demosaic::Bilinear, &IDENTITY); + let ks = kernels(phase, Demosaic::Bilinear, p.stride); + for y in 0..h { + for x in 0..w { + let i = (y + PAD) * p.stride + PAD + x; + for (c, k) in ks[(y & 1) * 2 + (x & 1)].iter().enumerate() { + let v: f32 = k.iter().map(|&(off, wt)| wt * p.data[(i as isize + off) as usize]).sum(); + assert!((fast[(y * w + x) * 3 + c] - v).abs() < 1e-6, "({x},{y}) c{c}"); + } + } + } + } + } + + #[test] + fn lab_table_is_accurate() { + let l = LabF::new(); + for i in 0..10000 { + let t = i as f32 / 4000.0; + assert!((l.f(t) - LabF::exact(t)).abs() < 2e-4, "{t}"); + } + } + + #[test] + fn tiny_planes_do_not_panic() { + for (w, h) in [(1, 1), (2, 1), (1, 3), (2, 2), (3, 3), (7, 9), (40, 1)] { + let cfa = vec![0.5f32; w * h]; + for m in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] { + assert_eq!(run(&cfa, w, h, [0, 1, 1, 2], m).len(), w * h * 3); + } + } + } +} diff --git a/crates/raw/src/develop.rs b/crates/raw/src/develop.rs new file mode 100644 index 00000000..dbfd8649 --- /dev/null +++ b/crates/raw/src/develop.rs @@ -0,0 +1,429 @@ +//! Raw development: linearize, subtract black, scale to white, white +//! balance, clip highlights, demosaic, camera → XYZ (D50) → linear ProPhoto, +//! exposure, gamma 1.8 encode to 16 bits, orientation. + +use crate::color::{self, Mat3}; +use crate::demosaic::{Demosaic, PAD, Padded, demosaic}; +use crate::error::{RawError, Result}; +use crate::sensor::Sensor; +use crate::{Limits, RawFormat, par}; + +/// White balance choice. +#[derive(Debug, Clone, Copy, PartialEq, Default)] +pub enum WhiteBalance { + /// The camera's as-shot white balance when the file records one, else automatic. + #[default] + AsShot, + /// Grey-world estimate from the image. + Auto, + /// Explicit R, G, B multipliers applied to the camera channels. + Multipliers([f64; 3]), +} + +/// Development settings. +#[derive(Debug, Clone, PartialEq)] +pub struct DevelopOptions { + pub demosaic: Demosaic, + pub white_balance: WhiteBalance, + /// Exposure compensation in EV, added to the file's baseline exposure. + pub exposure: f64, + /// Apply the file's orientation (rotate / flip the pixels). + pub orient: bool, + pub limits: Limits, +} + +impl Default for DevelopOptions { + fn default() -> Self { + DevelopOptions { demosaic: Demosaic::default(), white_balance: WhiteBalance::default(), exposure: 0.0, orient: true, limits: Limits::default() } + } +} + +/// What was developed and how. +#[derive(Debug, Clone, PartialEq)] +pub struct RawInfo { + pub format: RawFormat, + pub make: Option, + pub model: Option, + /// Sensor data size. + pub sensor_width: usize, + pub sensor_height: usize, + /// CFA pattern at the crop origin, e.g. `"RGGB"`; `None` for LinearRaw. + pub cfa: Option, + /// White-balance multipliers used (R, G, B, normalized to a minimum of 1). + pub wb_multipliers: [f64; 3], + pub orientation: u16, + pub baseline_exposure: f64, +} + +/// A developed image: 16-bit RGB in ProPhoto RGB (ROMM primaries, D50, gamma 1.8). +#[derive(Debug, Clone)] +pub struct Developed { + pub width: u32, + pub height: u32, + /// Interleaved RGB, `width * height * 3` samples. + pub rgb: Vec, + pub info: RawInfo, + /// Notes about anything approximated or not applied. + pub warnings: Vec, +} + +/// Linearization + black/white scaling for one sample (pre white balance), 0..1. +struct Scale<'a> { + s: &'a Sensor, + lin: Option<&'a [u16]>, +} + +impl Scale<'_> { + #[inline] + fn value(&self, x: usize, y: usize, sample: usize) -> f32 { + let s = self.s; + let raw = s.data[(y * s.width + x) * s.samples + sample]; + let v = match self.lin { + Some(l) => f32::from(l.get(usize::from(raw)).or_else(|| l.last()).copied().unwrap_or(raw)), + None => f32::from(raw), + }; + let black = s.black.at(x.saturating_sub(s.active.x), y.saturating_sub(s.active.y), sample, s.samples); + let white = s.white[sample.min(2)]; + let range = white - black; + if range.is_nan() || range < 1.0 { + return 0.0; + } + let n = (v - black) / range; + if n.is_finite() { n } else { 0.0 } + } +} + +/// Grey-world neutral (camera-space R, G, B) from unclipped mid-range samples. +fn grey_world(s: &Sensor, sc: &Scale) -> [f64; 3] { + let c = s.crop; + let step = (c.height / 512).max(1); + let rows: Vec = (0..c.height).step_by(step).collect(); + let parts = par::map(rows.len(), |i| { + let y = c.y + rows[i]; + let mut sum = [0.0f64; 3]; + let mut n = [0u64; 3]; + for x in c.x..c.x + c.width { + for k in 0..s.samples.min(3) { + let ch = match &s.cfa { + Some(cfa) => usize::from(cfa.color(x, y)), + None => k, + }; + let v = sc.value(x, y, k); + if v > 0.002 && v < 0.9 { + sum[ch] += f64::from(v); + n[ch] += 1; + } + } + } + (sum, n) + }); + let mut sum = [0.0f64; 3]; + let mut n = [0u64; 3]; + for (s, k) in parts { + for c in 0..3 { + sum[c] += s[c]; + n[c] += k[c]; + } + } + let mean = [0, 1, 2].map(|c| if n[c] > 0 { sum[c] / n[c] as f64 } else { 0.0 }); + if mean.iter().any(|m| m.is_nan() || *m <= 1e-6) { + return [1.0; 3]; + } + let m = mean.iter().cloned().fold(f64::MIN, f64::max); + mean.map(|v| v / m) +} + +fn normalize(n: [f64; 3]) -> Option<[f64; 3]> { + let m = n.iter().cloned().fold(f64::MIN, f64::max); + (m > 0.0 && n.iter().all(|v| v.is_finite() && *v > 1e-4)).then(|| n.map(|v| v / m)) +} + +/// Develops sensor data. +pub fn develop_sensor(s: &Sensor, opts: &DevelopOptions) -> Result { + let c = s.crop; + if c.is_empty() || c.x + c.width > s.width || c.y + c.height > s.height || s.data.len() < s.width * s.height * s.samples { + return Err(RawError::malformed("crop lies outside the sensor data")); + } + if s.samples != 1 && s.samples != 3 { + return Err(RawError::unsupported(format!("{} samples per pixel", s.samples))); + } + // Working buffers: one f32 plane plus f32 RGB plus u16 RGB per pixel. + opts.limits.check(c.width as u64, c.height as u64, 4 + 12 + 6)?; + let mut warnings = s.warnings.clone(); + let sc = Scale { s, lin: s.linearization.as_deref() }; + + // White balance: camera-space neutral, max component 1. + let as_shot = s + .color + .as_shot_neutral + .or_else(|| s.color.as_shot_white_xy.and_then(|xy| s.color.xy_to_neutral(xy))) + .or_else(|| s.camera_wb.map(|m| m.map(|v| 1.0 / v))) + .and_then(normalize); + let neutral = match opts.white_balance { + WhiteBalance::AsShot => match as_shot { + Some(n) => n, + None => { + warnings.push("no as-shot white balance in the file; estimated automatically (grey world)".into()); + grey_world(s, &sc) + } + }, + WhiteBalance::Auto => grey_world(s, &sc), + WhiteBalance::Multipliers(m) => { + normalize(m.map(|v| 1.0 / v)).ok_or_else(|| RawError::Malformed("white-balance multipliers must be positive".into()))? + } + }; + let mult = neutral.map(|v| (1.0 / v) as f32); + + // Balanced camera → XYZ D50. + let to_xyz: Mat3 = match s.color.balanced_to_xyz_d50(neutral) { + Some(m) => m, + None => { + if s.format != RawFormat::Dng { + warnings.push("no colour calibration for this camera; camera colours are treated as sRGB primaries".into()); + } + color::srgb_to_xyz_d50() + } + }; + let gain = 2f64.powf(s.baseline_exposure + opts.exposure.clamp(-10.0, 10.0)); + let out_m = color::scale(&color::mul(&color::xyz_d50_to_prophoto(), &to_xyz), gain); + let out_m: [[f32; 3]; 3] = out_m.map(|r| r.map(|v| v as f32)); + + let (w, h) = (c.width, c.height); + let band = par::band_rows(w); + let mut rgb: Vec; + let phase; + match &s.cfa { + Some(cfa) => { + // Normalized, balanced, clipped CFA plane. + let mut plane = Padded::new(w, h); + let stride = plane.stride; + // Fast path: black varies at most per 2×2 position (the Bayer period). + let bl = &s.black; + let periodic = bl.delta_h.is_empty() && bl.delta_v.is_empty() && matches!(bl.rows, 1 | 2) && matches!(bl.cols, 1 | 2); + plane.fill_rows(band, |y0, chunk| { + for (r, prow) in chunk.chunks_exact_mut(stride).enumerate() { + let row = &mut prow[PAD..PAD + w]; + let y = c.y + y0 + r; + if periodic { + // Per column parity: black and the combined scale × white-balance gain. + let k: [(f32, f32); 2] = std::array::from_fn(|p| { + let x = c.x + p; + let black = bl.at(x.saturating_sub(s.active.x), y.saturating_sub(s.active.y), 0, 1); + let range = s.white[0] - black; + let gain = if range >= 1.0 { mult[usize::from(cfa.color(x, y)).min(2)] / range } else { 0.0 }; + (black, gain) + }); + let src = &s.data[y * s.width + c.x..y * s.width + c.x + w]; + match &sc.lin { + None => { + for (i, (v, &raw)) in row.iter_mut().zip(src).enumerate() { + let (black, gain) = k[i & 1]; + *v = ((f32::from(raw) - black) * gain).clamp(0.0, 1.0); + } + } + Some(l) => { + for (i, (v, &raw)) in row.iter_mut().zip(src).enumerate() { + let (black, gain) = k[i & 1]; + let lin = f32::from(l.get(usize::from(raw)).or_else(|| l.last()).copied().unwrap_or(raw)); + *v = ((lin - black) * gain).clamp(0.0, 1.0); + } + } + } + } else { + for (i, v) in row.iter_mut().enumerate() { + let x = c.x + i; + let ch = usize::from(cfa.color(x, y)).min(2); + *v = (sc.value(x, y, 0) * mult[ch]).clamp(0.0, 1.0); + } + } + // Lens-shading gain maps, after the white-balance clip so clipped + // highlights stay neutral. + if let Some(g) = gain_row(s, y, c.x, w, 0) { + for (v, g) in row.iter_mut().zip(g) { + *v = (*v * g).clamp(0.0, 1.0); + } + } + } + }); + plane.fill_borders(); + let p = cfa.phase(c.x, c.y); + phase = Some(p); + rgb = demosaic(&plane, p, opts.demosaic, &to_xyz); + drop(plane); + } + None => { + phase = None; + rgb = vec![0.0f32; w * h * 3]; + par::chunks_mut(&mut rgb, band * w * 3, |b, chunk| { + for (r, row) in chunk.chunks_exact_mut(w * 3).enumerate() { + let y = c.y + b * band + r; + for i in 0..w { + for k in 0..3 { + row[i * 3 + k] = (sc.value(c.x + i, y, k) * mult[k]).clamp(0.0, 1.0); + } + } + for k in 0..3 { + if let Some(g) = gain_row(s, y, c.x, w, k) { + for (i, g) in g.into_iter().enumerate() { + row[i * 3 + k] = (row[i * 3 + k] * g).clamp(0.0, 1.0); + } + } + } + } + }); + } + } + + // Camera → ProPhoto, exposure, gamma 1.8, 16 bits. + let mut out = vec![0u16; w * h * 3]; + let lut = GammaLut::new(); + par::chunks_mut(&mut out, band * w * 3, |b, chunk| { + let start = b * band * w * 3; + let src = &rgb[start..start + chunk.len()]; + for (o, p) in chunk.as_chunks_mut::<3>().0.iter_mut().zip(src.as_chunks::<3>().0) { + for k in 0..3 { + o[k] = lut.encode(out_m[k][0] * p[0] + out_m[k][1] * p[1] + out_m[k][2] * p[2]); + } + } + }); + rgb.clear(); + rgb.shrink_to_fit(); + + let (out, ow, oh) = if opts.orient { orient(out, w, h, s.orientation) } else { (out, w, h) }; + let nmax = mult.iter().cloned().fold(f32::MAX, f32::min).max(1e-6); + Ok(Developed { + width: ow as u32, + height: oh as u32, + rgb: out, + info: RawInfo { + format: s.format, + make: s.make.clone(), + model: s.model.clone(), + sensor_width: s.width, + sensor_height: s.height, + cfa: phase.map(|p| p.iter().map(|c| ['R', 'G', 'B'][usize::from(*c).min(2)]).collect()), + wb_multipliers: mult.map(|m| f64::from(m / nmax)), + orientation: s.orientation, + baseline_exposure: s.baseline_exposure, + }, + warnings, + }) +} + +/// Gain-map factors (DNG OpcodeList2 GainMap) for data row `y`, columns +/// `x0..x0 + w`, sample `k`; `None` when no map touches the row. +fn gain_row(s: &Sensor, y: usize, x0: usize, w: usize, k: usize) -> Option> { + let a = s.active; + if s.gain_maps.is_empty() || y < a.y || a.width == 0 || a.height == 0 { + return None; + } + let ay = y - a.y; + // Map points are placed in coordinates relative to the (active) image. + let rv = (ay as f64 + 0.5) / a.height as f64; + let mut out: Option> = None; + for m in &s.gain_maps { + if k < m.plane || k >= m.plane.saturating_add(m.planes) || ay < m.top || ay >= m.bottom || !(ay - m.top).is_multiple_of(m.row_pitch) { + continue; + } + let row = m.row(rv, k - m.plane); + let o = out.get_or_insert_with(|| vec![1.0; w]); + for (i, f) in o.iter_mut().enumerate() { + let x = x0 + i; + if x >= a.x && m.covers(x - a.x, ay) { + *f *= m.at(&row, ((x - a.x) as f64 + 0.5) / a.width as f64); + } + } + } + out +} + +/// Linear 0..1 → gamma-1.8 16-bit code. Tabulated over `s = sqrt(v)`, where +/// the curve `s^(2/1.8)` is smooth down to 0, with linear interpolation +/// (error below 0.2 code values). +struct GammaLut { + table: Vec, +} + +const GAMMA_STEPS: usize = 4096; + +impl GammaLut { + fn new() -> Self { + GammaLut { table: (0..=GAMMA_STEPS + 1).map(|i| ((i as f32 / GAMMA_STEPS as f32).min(1.0)).powf(2.0 / 1.8) * 65535.0).collect() } + } + + #[inline] + fn encode(&self, v: f32) -> u16 { + let v = if v.is_nan() { 0.0 } else { v.clamp(0.0, 1.0) }; + let x = v.sqrt() * GAMMA_STEPS as f32; + let i = (x as usize).min(GAMMA_STEPS); + let f = x - i as f32; + let (a, b) = (self.table[i], self.table[i + 1]); + (a + (b - a) * f + 0.5) as u16 + } +} + +/// Applies a TIFF orientation (1–8) to interleaved RGB. +pub(crate) fn orient(src: Vec, w: usize, h: usize, o: u16) -> (Vec, usize, usize) { + if !(2..=8).contains(&o) { + return (src, w, h); + } + let swap = o >= 5; + let (ow, oh) = if swap { (h, w) } else { (w, h) }; + let mut out = vec![0u16; src.len()]; + let band = par::band_rows(ow); + par::chunks_mut(&mut out, band * ow * 3, |b, chunk| { + for (r, row) in chunk.chunks_exact_mut(ow * 3).enumerate() { + let y2 = b * band + r; + for x2 in 0..ow { + let (x, y) = match o { + 2 => (w - 1 - x2, y2), + 3 => (w - 1 - x2, h - 1 - y2), + 4 => (x2, h - 1 - y2), + 5 => (y2, x2), + 6 => (y2, h - 1 - x2), + 7 => (w - 1 - y2, h - 1 - x2), + _ => (w - 1 - y2, x2), + }; + let i = (y * w + x) * 3; + row[x2 * 3..x2 * 3 + 3].copy_from_slice(&src[i..i + 3]); + } + } + }); + (out, ow, oh) +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn gamma_table_matches_powf() { + let lut = GammaLut::new(); + for i in 0..=100_000 { + let v = i as f32 / 100_000.0; + let exact = v.powf(1.0 / 1.8) * 65535.0; + assert!((f32::from(lut.encode(v)) - exact).abs() <= 1.0, "{v}: {} vs {exact}", lut.encode(v)); + } + assert_eq!(lut.encode(-1.0), 0); + assert_eq!(lut.encode(2.0), 65535); + assert_eq!(lut.encode(f32::NAN), 0); + } + + #[test] + fn orientations() { + // 3×2 image, pixel value = index. + let src: Vec = (0..6).flat_map(|i| [i, i, i]).collect(); + let px = |v: &[u16]| v.as_chunks::<3>().0.iter().map(|p| p[0]).collect::>(); + assert_eq!(px(&orient(src.clone(), 3, 2, 1).0), vec![0, 1, 2, 3, 4, 5]); + assert_eq!(px(&orient(src.clone(), 3, 2, 2).0), vec![2, 1, 0, 5, 4, 3]); + assert_eq!(px(&orient(src.clone(), 3, 2, 3).0), vec![5, 4, 3, 2, 1, 0]); + assert_eq!(px(&orient(src.clone(), 3, 2, 4).0), vec![3, 4, 5, 0, 1, 2]); + let (v, w, h) = orient(src.clone(), 3, 2, 6); + assert_eq!((w, h), (2, 3)); + // Rotated 90° clockwise: the left column (0, 3) becomes the top row (3, 0). + assert_eq!(px(&v), vec![3, 0, 4, 1, 5, 2]); + assert_eq!(px(&orient(src.clone(), 3, 2, 8).0), vec![2, 5, 1, 4, 0, 3]); + assert_eq!(px(&orient(src.clone(), 3, 2, 5).0), vec![0, 3, 1, 4, 2, 5]); + assert_eq!(px(&orient(src, 3, 2, 7).0), vec![5, 2, 4, 1, 3, 0]); + } +} diff --git a/crates/raw/src/dng.rs b/crates/raw/src/dng.rs new file mode 100644 index 00000000..6ae428cd --- /dev/null +++ b/crates/raw/src/dng.rs @@ -0,0 +1,225 @@ +//! Digital Negative (Adobe DNG Specification 1.7): raw IFD selection, CFA +//! description, linearization, black/white levels, active area, default +//! crop, colour calibration and white balance. + +use crate::color::{self, Calibration, ColorInfo, IDENTITY}; +use crate::error::{RawError, Result}; +use crate::sensor::{BlackLevels, Cfa, JpegLayout, Rect, Sensor, read_plane}; +use crate::tiff::{Ifd, Tiff, tag}; +use crate::{Limits, RawFormat, opcodes}; + +const PHOTOMETRIC_CFA: u32 = 32803; +const PHOTOMETRIC_LINEAR_RAW: u32 = 34892; + +/// The main raw image: NewSubFileType 0 with CFA or LinearRaw data, largest first. +fn raw_ifd(t: &Tiff, ifds: &[Ifd]) -> Option { + ifds.iter() + .filter(|i| t.tag_uint(i, tag::NEW_SUBFILE_TYPE).unwrap_or(0) == 0) + .filter(|i| matches!(t.tag_uint(i, tag::PHOTOMETRIC), Some(PHOTOMETRIC_CFA | PHOTOMETRIC_LINEAR_RAW))) + .max_by_key(|i| u64::from(t.tag_uint(i, tag::IMAGE_WIDTH).unwrap_or(0)) * u64::from(t.tag_uint(i, tag::IMAGE_LENGTH).unwrap_or(0))) + .cloned() +} + +fn matrix(t: &Tiff, ifd: &Ifd, tg: u16) -> Option { + let v = t.tag_floats(ifd, tg); + if v.len() != 9 { + return None; + } + color::from_slice(&v) +} + +fn triple(t: &Tiff, ifd: &Ifd, tg: u16) -> Option<[f64; 3]> { + let v = t.tag_floats(ifd, tg); + (v.len() == 3 && v.iter().all(|x| x.is_finite() && *x > 0.0)).then(|| [v[0], v[1], v[2]]) +} + +/// Colour calibration from IFD0. +pub(crate) fn color_info(t: &Tiff, ifd0: &Ifd) -> ColorInfo { + let mut cal = Vec::new(); + for (cm, cc, fm, ill, default_t) in [ + (tag::COLOR_MATRIX_1, tag::CAMERA_CALIBRATION_1, tag::FORWARD_MATRIX_1, tag::CALIBRATION_ILLUMINANT_1, 2856.0), + (tag::COLOR_MATRIX_2, tag::CAMERA_CALIBRATION_2, tag::FORWARD_MATRIX_2, tag::CALIBRATION_ILLUMINANT_2, 6504.0), + ] { + let Some(color_matrix) = matrix(t, ifd0, cm) else { continue }; + if color::invert(&color_matrix).is_none() { + continue; + } + let temperature = t.tag_uint(ifd0, ill).and_then(color::illuminant_temperature).unwrap_or(default_t); + cal.push(Calibration { + temperature, + color_matrix, + forward_matrix: matrix(t, ifd0, fm).filter(|m| color::invert(m).is_some()), + camera_calibration: matrix(t, ifd0, cc).filter(|m| color::invert(m).is_some()).unwrap_or(IDENTITY), + }); + } + cal.sort_by(|a, b| a.temperature.total_cmp(&b.temperature)); + if cal.len() == 2 && (cal[0].forward_matrix.is_some() != cal[1].forward_matrix.is_some()) { + for c in &mut cal { + c.forward_matrix = None; + } + } + let white_xy = { + let v = t.tag_floats(ifd0, tag::AS_SHOT_WHITE_XY); + (v.len() == 2 && v.iter().all(|x| *x > 0.0 && *x < 1.0)).then(|| [v[0], v[1]]) + }; + ColorInfo { + calibrations: cal, + analog_balance: triple(t, ifd0, tag::ANALOG_BALANCE), + as_shot_neutral: triple(t, ifd0, tag::AS_SHOT_NEUTRAL), + as_shot_white_xy: white_xy, + } +} + +/// CFA pattern of a DNG raw IFD, anchored at the active area origin. +fn cfa(t: &Tiff, ifd: &Ifd, active: &Rect) -> Result { + let dim = t.tag_uints(ifd, tag::CFA_REPEAT_PATTERN_DIM); + let (rows, cols) = match dim.as_slice() { + [r, c] => (*r as usize, *c as usize), + _ => (2, 2), + }; + let pattern = t.tag_uints(ifd, tag::CFA_PATTERN); + if rows == 0 || cols == 0 || rows > 16 || cols > 16 || pattern.len() != rows * cols { + return Err(RawError::malformed("bad CFA pattern")); + } + let planes = match t.tag_uints(ifd, tag::CFA_PLANE_COLOR) { + v if v.is_empty() => vec![0, 1, 2], + v => v, + }; + if planes != [0, 1, 2] { + return Err(RawError::unsupported(format!("CFA plane colours {planes:?} (only red, green, blue)"))); + } + if t.tag_uint(ifd, tag::CFA_LAYOUT).unwrap_or(1) != 1 { + return Err(RawError::unsupported("non-rectangular CFA layout")); + } + let colors = pattern + .into_iter() + .map(|code| u8::try_from(code).ok().filter(|c| *c <= 2)) + .collect::>>() + .ok_or_else(|| RawError::unsupported("CFA colours other than red, green and blue"))?; + let c = Cfa { width: cols, height: rows, colors, origin_x: active.x, origin_y: active.y }; + if !c.is_bayer() { + return Err(RawError::unsupported(format!("{rows}x{cols} non-Bayer CFA (e.g. X-Trans)"))); + } + Ok(c) +} + +fn rect_from(t: &Tiff, ifd: &Ifd, width: usize, height: usize) -> Rect { + let full = Rect::new(0, 0, width, height); + match t.tag_uints(ifd, tag::ACTIVE_AREA).as_slice() { + [top, left, bottom, right] if bottom > top && right > left => { + Rect::new(*left as usize, *top as usize, (*right - *left) as usize, (*bottom - *top) as usize).intersect(&full) + } + _ => full, + } +} + +/// Decodes a DNG file's raw image and metadata. +pub(crate) fn decode(t: &Tiff, limits: &Limits) -> Result { + let ifds = t.all_ifds(); + let ifd0 = t.ifd_at(t.first_ifd, 0).ok_or_else(|| RawError::malformed("DNG has no IFD0"))?; + let raw = raw_ifd(t, &ifds).ok_or_else(|| RawError::malformed("DNG has no raw image"))?; + let mut warnings = Vec::new(); + let photometric = t.tag_uint(&raw, tag::PHOTOMETRIC).unwrap_or(0); + let plane = read_plane(t, &raw, limits, JpegLayout::Flat)?; + let samples = plane.samples; + let is_cfa = photometric == PHOTOMETRIC_CFA; + if is_cfa && samples != 1 { + return Err(RawError::unsupported(format!("CFA data with {samples} samples per pixel"))); + } + if !is_cfa && samples != 3 { + return Err(RawError::unsupported(format!("LinearRaw data with {samples} samples per pixel"))); + } + let active = rect_from(t, &raw, plane.width, plane.height); + if active.is_empty() { + return Err(RawError::malformed("empty active area")); + } + let cfa = if is_cfa { Some(cfa(t, &raw, &active)?) } else { None }; + + // Default crop, relative to the active area. + let mut crop = active; + let origin = t.tag_floats(&raw, tag::DEFAULT_CROP_ORIGIN); + let size = t.tag_floats(&raw, tag::DEFAULT_CROP_SIZE); + if let ([ox, oy], [cw, ch]) = (origin.as_slice(), size.as_slice()) { + let f = |v: f64| if v.is_finite() && v >= 0.0 { v.round() as usize } else { 0 }; + let r = Rect::new(active.x + f(*ox), active.y + f(*oy), f(*cw), f(*ch)).intersect(&active); + if !r.is_empty() { + crop = r; + } + } + let scale = t.tag_floats(&raw, tag::DEFAULT_SCALE); + if scale.len() == 2 && (scale[0] - scale[1]).abs() > 1e-6 { + warnings.push("non-square pixels (DefaultScale) are not resampled".to_string()); + } + // Opcode lists: OpcodeList2 gain maps (lens shading) are applied; everything else is reported. + let mut gain_maps = Vec::new(); + for (tg, what) in [(tag::OPCODE_LIST_1, 1), (tag::OPCODE_LIST_2, 2), (tag::OPCODE_LIST_3, 3)] { + let Some(bytes) = raw.get(tg).and_then(|e| t.raw(e)) else { continue }; + let (maps, skipped) = opcodes::parse(bytes); + if what == 2 { + gain_maps = maps; + } else if !maps.is_empty() { + warnings.push(format!("DNG OpcodeList{what}: GainMap is not applied")); + } + if !skipped.is_empty() { + warnings.push(format!("DNG OpcodeList{what}: {} not applied", skipped.join(", "))); + } + } + + let max = ((1u32 << plane.bits) - 1) as f32; + let linearization = raw.get(tag::LINEARIZATION_TABLE).map(|e| t.uints(e)).filter(|v| !v.is_empty() && v.len() <= 65536); + let linearization: Option> = linearization.map(|v| v.into_iter().map(|x| x.min(65535) as u16).collect()); + let white_default = match &linearization { + Some(l) => f32::from(l.iter().copied().max().unwrap_or(u16::MAX)), + None => max, + }; + let wl = t.tag_floats(&raw, tag::WHITE_LEVEL); + let white = [0, 1, 2].map(|s| wl.get(s).or_else(|| wl.first()).map(|v| *v as f32).filter(|v| *v > 0.0).unwrap_or(white_default)); + + let (rows, cols) = match t.tag_uints(&raw, tag::BLACK_LEVEL_REPEAT_DIM).as_slice() { + [r, c] if (1..=16).contains(r) && (1..=16).contains(c) => (*r as usize, *c as usize), + _ => (1, 1), + }; + let bl = t.tag_floats(&raw, tag::BLACK_LEVEL); + let mut black = if bl.len() == rows * cols * samples { + BlackLevels { rows, cols, values: bl.iter().map(|v| *v as f32).collect(), delta_h: Vec::new(), delta_v: Vec::new() } + } else { + BlackLevels::uniform(bl.first().copied().unwrap_or(0.0) as f32) + }; + let dh = t.tag_floats(&raw, tag::BLACK_LEVEL_DELTA_H); + if dh.len() == active.width { + black.delta_h = dh.iter().map(|v| *v as f32).collect(); + } + let dv = t.tag_floats(&raw, tag::BLACK_LEVEL_DELTA_V); + if dv.len() == active.height { + black.delta_v = dv.iter().map(|v| *v as f32).collect(); + } + + let color = color_info(t, &ifd0); + if color.calibrations.is_empty() { + warnings.push("the DNG has no colour matrix; colours are approximate".to_string()); + } + let baseline_exposure = t.tag_floats(&ifd0, tag::BASELINE_EXPOSURE).first().copied().unwrap_or(0.0) + + t.tag_floats(&raw, tag::BASELINE_EXPOSURE_OFFSET).first().copied().unwrap_or(0.0); + let model = t.tag_ascii(&ifd0, tag::UNIQUE_CAMERA_MODEL).or_else(|| t.tag_ascii(&ifd0, tag::MODEL)); + Ok(Sensor { + format: RawFormat::Dng, + make: t.tag_ascii(&ifd0, tag::MAKE), + model, + width: plane.width, + height: plane.height, + samples, + data: plane.data, + cfa, + linearization, + black, + white, + active, + crop, + color, + camera_wb: None, + orientation: t.tag_uint(&ifd0, tag::ORIENTATION).map(|o| o as u16).filter(|o| (1..=8).contains(o)).unwrap_or(1), + baseline_exposure: if baseline_exposure.is_finite() { baseline_exposure.clamp(-10.0, 10.0) } else { 0.0 }, + gain_maps, + warnings, + }) +} diff --git a/crates/raw/src/error.rs b/crates/raw/src/error.rs new file mode 100644 index 00000000..ff9a1dc1 --- /dev/null +++ b/crates/raw/src/error.rs @@ -0,0 +1,40 @@ +//! Error type. + +use std::fmt; + +/// Why a raw file could not be decoded. +#[derive(Debug, Clone, PartialEq, Eq)] +pub enum RawError { + /// Not a camera raw file this crate recognizes. + NotRaw, + /// A raw file (or a variant of one) that is recognized but not decoded yet. + Unsupported(String), + /// The file is damaged or inconsistent. + Malformed(String), + /// Declared sizes exceed the [`crate::Limits`]. + LimitExceeded(String), +} + +impl RawError { + pub(crate) fn malformed(msg: impl Into) -> Self { + RawError::Malformed(msg.into()) + } + pub(crate) fn unsupported(msg: impl Into) -> Self { + RawError::Unsupported(msg.into()) + } +} + +impl fmt::Display for RawError { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + match self { + RawError::NotRaw => write!(f, "not a camera raw file"), + RawError::Unsupported(m) => write!(f, "unsupported camera raw: {m}"), + RawError::Malformed(m) => write!(f, "malformed camera raw: {m}"), + RawError::LimitExceeded(m) => write!(f, "camera raw exceeds limits: {m}"), + } + } +} + +impl std::error::Error for RawError {} + +pub(crate) type Result = std::result::Result; diff --git a/crates/raw/src/lib.rs b/crates/raw/src/lib.rs new file mode 100644 index 00000000..8871e319 --- /dev/null +++ b/crates/raw/src/lib.rs @@ -0,0 +1,234 @@ +//! `photocraft-raw`: a clean-room, pure-Rust camera raw decoder and developer. +//! +//! Implemented only from public specifications and papers: TIFF 6.0, +//! TIFF/EP (ISO 12234-2), the Adobe DNG Specification 1.7, ITU-T T.81 +//! (lossless JPEG, process 14), the published structure of Canon's CR2 +//! container, and the demosaicing papers cited in [`Demosaic`]. No code from +//! dcraw, LibRaw, rawspeed, rawler, rawloader or darktable was used. +//! +//! * [`identify`] / [`is_raw`] recognise raw files from their bytes. +//! * [`decode`] reads the undeveloped [`Sensor`] data and metadata. +//! * [`develop`] turns it into 16-bit RGB in ProPhoto RGB (ROMM primaries, +//! D50 white, gamma 1.8): linearization, black/white levels, white balance, +//! highlight clipping, demosaicing, camera → XYZ → ProPhoto, exposure and +//! orientation. +//! * [`embedded_preview`] finds the camera's full-size JPEG preview, for +//! formats whose sensor data is not decoded yet. +//! +//! Decoded today: DNG (uncompressed and lossless-JPEG, strips and tiles, CFA +//! and LinearRaw), CR2 (lossless JPEG with Canon slices) and uncompressed or +//! lossless-JPEG TIFF/EP raws (NEF, ARW, PEF… when not vendor-compressed). +//! Everything else reports [`RawError::Unsupported`]. +//! +//! The crate is standalone (no workspace dependencies), does no I/O, builds for +//! `wasm32-unknown-unknown` and never panics on hostile input: sizes are +//! checked against [`Limits`] before allocating. + +#![forbid(unsafe_code)] +#![deny(clippy::unwrap_used, clippy::expect_used, clippy::panic, clippy::unimplemented, clippy::todo, clippy::unreachable)] + +mod color; +mod cr2; +mod demosaic; +mod develop; +mod dng; +mod error; +mod ljpeg; +mod opcodes; +mod par; +mod preview; +mod sensor; +mod tiff; +mod tiffep; + +#[cfg(feature = "testgen")] +pub mod testgen; + +pub use crate::color::{Calibration, ColorInfo, Mat3}; +pub use crate::demosaic::Demosaic; +pub use crate::develop::{DevelopOptions, Developed, RawInfo, WhiteBalance, develop_sensor}; +pub use crate::error::RawError; +pub use crate::opcodes::GainMap; +pub use crate::preview::{Preview, embedded_preview}; +pub use crate::sensor::{BlackLevels, Cfa, Rect, Sensor}; + +use crate::tiff::{Tiff, tag}; + +/// Name of the output colour space of [`develop`]. +pub const OUTPUT_SPACE: &str = "ProPhoto RGB (ROMM primaries, D50, gamma 1.8)"; + +/// Camera raw container families. +#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)] +pub enum RawFormat { + Dng, + /// Canon CR2. + Cr2, + /// Canon CR3 (ISO BMFF). + Cr3, + /// Nikon NEF / NRW. + Nef, + /// Sony ARW / SR2. + Arw, + /// Pentax PEF. + Pef, + /// Olympus ORF. + Orf, + /// Panasonic RW2. + Rw2, + /// Fujifilm RAF. + Raf, + /// Another TIFF/EP raw. + TiffEp, +} + +impl RawFormat { + pub fn name(self) -> &'static str { + match self { + RawFormat::Dng => "DNG", + RawFormat::Cr2 => "CR2", + RawFormat::Cr3 => "CR3", + RawFormat::Nef => "NEF", + RawFormat::Arw => "ARW", + RawFormat::Pef => "PEF", + RawFormat::Orf => "ORF", + RawFormat::Rw2 => "RW2", + RawFormat::Raf => "RAF", + RawFormat::TiffEp => "TIFF/EP raw", + } + } +} + +/// File extensions of the raw formats this crate recognises. +pub const EXTENSIONS: &[&str] = &["dng", "cr2", "cr3", "nef", "nrw", "arw", "srf", "sr2", "pef", "orf", "rw2", "raf"]; + +/// Decompression-bomb guards, checked on header values before allocating +/// (the same shape and defaults as `photocraft-codecs`' limits). +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub struct Limits { + pub max_width: u32, + pub max_height: u32, + /// Maximum `width * height`. + pub max_pixels: u64, + /// Maximum bytes of one buffer. + pub max_alloc: u64, +} + +impl Default for Limits { + fn default() -> Self { + Limits { max_width: 1 << 18, max_height: 1 << 18, max_pixels: 1 << 28, max_alloc: 2 << 30 } + } +} + +impl Limits { + /// Checks a `width × height` buffer of `bytes_per_pixel` bytes per pixel. + pub fn check(&self, width: u64, height: u64, bytes_per_pixel: u64) -> Result<(), RawError> { + if width == 0 || height == 0 { + return Err(RawError::Malformed(format!("zero-sized image {width}x{height}"))); + } + if width > u64::from(self.max_width) || height > u64::from(self.max_height) { + return Err(RawError::LimitExceeded(format!("dimensions {width}x{height} exceed {}x{}", self.max_width, self.max_height))); + } + let pixels = width.saturating_mul(height); + if pixels > self.max_pixels { + return Err(RawError::LimitExceeded(format!("{pixels} pixels exceed {}", self.max_pixels))); + } + let bytes = pixels.saturating_mul(bytes_per_pixel); + if bytes > self.max_alloc || usize::try_from(bytes).is_err() { + return Err(RawError::LimitExceeded(format!("{bytes} bytes exceed {}", self.max_alloc))); + } + Ok(()) + } +} + +const PHOTOMETRIC_CFA: u32 = 32803; +const PHOTOMETRIC_LINEAR_RAW: u32 = 34892; + +fn make_format(make: Option<&str>) -> RawFormat { + let m = make.unwrap_or("").to_ascii_uppercase(); + if m.starts_with("NIKON") { + RawFormat::Nef + } else if m.starts_with("SONY") { + RawFormat::Arw + } else if m.starts_with("PENTAX") || m.starts_with("RICOH") { + RawFormat::Pef + } else if m.starts_with("OLYMPUS") || m.starts_with("OM DIGITAL") { + RawFormat::Orf + } else { + RawFormat::TiffEp + } +} + +/// Recognises a camera raw file from its bytes (cheap: header and IFDs only). +pub fn identify(bytes: &[u8]) -> Option { + if bytes.get(4..12) == Some(b"ftypcrx ") { + return Some(RawFormat::Cr3); + } + if bytes.starts_with(b"FUJIFILMCCD-RAW") { + return Some(RawFormat::Raf); + } + let t = Tiff::new(bytes)?; + match bytes.get(0..4) { + Some(b"IIRO" | b"IIRS" | b"MMOR") => return Some(RawFormat::Orf), + Some(b"IIU\0") => return Some(RawFormat::Rw2), + _ => {} + } + if bytes.get(8..11) == Some(b"CR\x02") { + return Some(RawFormat::Cr2); + } + let ifds = t.all_ifds(); + if ifds.iter().any(|i| i.has(tag::DNG_VERSION)) { + return Some(RawFormat::Dng); + } + if ifds.iter().any(|i| matches!(t.tag_uint(i, tag::PHOTOMETRIC), Some(PHOTOMETRIC_CFA | PHOTOMETRIC_LINEAR_RAW))) { + let make = ifds.iter().find_map(|i| t.tag_ascii(i, tag::MAKE)); + return Some(make_format(make.as_deref())); + } + if tiffep::is_sony_non_cfa(&t, &ifds) { + return Some(RawFormat::Arw); + } + None +} + +/// `true` if `bytes` look like a camera raw file. +pub fn is_raw(bytes: &[u8]) -> bool { + identify(bytes).is_some() +} + +/// Reads the undeveloped sensor data and metadata. +pub fn decode(bytes: &[u8], limits: &Limits) -> Result { + let format = identify(bytes).ok_or(RawError::NotRaw)?; + match format { + RawFormat::Cr3 => Err(RawError::unsupported("Canon CR3 (ISO BMFF / CRX) is not decoded yet")), + RawFormat::Raf => Err(RawError::unsupported("Fujifilm RAF is not decoded yet")), + RawFormat::Rw2 => Err(RawError::unsupported("Panasonic RW2 is not decoded yet")), + RawFormat::Orf => Err(RawError::unsupported("Olympus ORF is not decoded yet")), + _ => { + let t = Tiff::new(bytes).ok_or(RawError::NotRaw)?; + match format { + RawFormat::Dng => dng::decode(&t, limits), + RawFormat::Cr2 => cr2::decode(&t, limits), + f => tiffep::decode(&t, f, limits), + } + } + } +} + +/// Decodes a lossless JPEG (T.81 process 14) stream into `(width, height, +/// components, samples)`; at most `max_samples` samples are allocated. +/// Exposed for fuzzing and tools. +#[doc(hidden)] +pub fn decode_lossless_jpeg(bytes: &[u8], max_samples: usize) -> Result<(usize, usize, usize, Vec), RawError> { + ljpeg::decode(bytes, max_samples).map(|(f, v)| (f.width, f.height, f.components, v)) +} + +/// Lists the TIFF structure of a raw file (debugging aid). +#[doc(hidden)] +pub fn dump_structure(bytes: &[u8]) -> String { + tiff::dump(bytes) +} + +/// Decodes and develops a raw file. +pub fn develop(bytes: &[u8], opts: &DevelopOptions) -> Result { + let sensor = decode(bytes, &opts.limits)?; + develop_sensor(&sensor, opts) +} diff --git a/crates/raw/src/ljpeg.rs b/crates/raw/src/ljpeg.rs new file mode 100644 index 00000000..e9c42bde --- /dev/null +++ b/crates/raw/src/ljpeg.rs @@ -0,0 +1,576 @@ +//! Lossless JPEG decoder (ITU-T T.81 / ISO 10918-1, process 14 "LJ92": +//! Huffman-coded, predictive, sequential, 2–16 bit samples). +//! +//! Written from the T.81 specification: Annex C (Huffman table +//! specification), Annex F.2.2.1-3 (decoding of DC-style difference +//! categories, here per Annex H) and Annex H (lossless mode: predictors 1–7, +//! the point transform, first-line / first-column rules, restart intervals). +//! +//! Output is the sample sequence in raster order: `height` rows of `width × +//! components` interleaved samples. Only one interleaved scan holding every +//! component with sampling factors 1×1 is supported, which is how DNG, CR2 +//! and the TIFF-based raws use it. + +use crate::error::{RawError, Result}; + +/// Frame header (SOF3). +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub(crate) struct Frame { + pub precision: u8, + pub width: usize, + pub height: usize, + pub components: usize, +} + +impl Frame { + pub fn samples(&self) -> Option { + self.width.checked_mul(self.height)?.checked_mul(self.components) + } +} + +/// Bits resolved by one table lookup. +const LOOKUP_BITS: u32 = 13; + +#[derive(Clone)] +struct Huffman { + /// `(length << 8) | symbol` for codes up to `LOOKUP_BITS` long, 0 if longer. + lookup: Vec, + /// `(difference << 8) | bits` when a code and its extra bits fit in + /// `LOOKUP_BITS`, else 0: one lookup per sample on the common path. + fast: Vec, + /// Largest code of each length (-1 = none), lengths 1..=16. + maxcode: [i32; 17], + mincode: [i32; 17], + valptr: [usize; 17], + vals: Vec, +} + +impl Huffman { + /// Builds the decoder from BITS (counts per length 1..16) and HUFFVAL (T.81 Annex C). + fn new(counts: &[u8; 16], vals: &[u8]) -> Result { + let total: usize = counts.iter().map(|&c| usize::from(c)).sum(); + if total == 0 || total > 256 || vals.len() < total { + return Err(RawError::malformed("bad Huffman table")); + } + let mut h = Huffman { + lookup: vec![0; 1 << LOOKUP_BITS], + fast: vec![0; 1 << LOOKUP_BITS], + maxcode: [-1; 17], + mincode: [0; 17], + valptr: [0; 17], + vals: vals[..total].to_vec(), + }; + let mut code: u32 = 0; + let mut k = 0usize; + for len in 1..=16u32 { + let n = usize::from(counts[len as usize - 1]); + if n > 0 { + h.valptr[len as usize] = k; + h.mincode[len as usize] = code as i32; + for _ in 0..n { + if code >= (1 << len) { + return Err(RawError::malformed("Huffman code overflow")); + } + if len <= LOOKUP_BITS { + let shift = LOOKUP_BITS - len; + let base = (code << shift) as usize; + let entry = ((len as u16) << 8) | u16::from(h.vals[k]); + for slot in &mut h.lookup[base..base + (1 << shift)] { + *slot = entry; + } + } + code += 1; + k += 1; + } + h.maxcode[len as usize] = code as i32 - 1; + } + code <<= 1; + } + for (idx, f) in h.fast.iter_mut().enumerate() { + let e = h.lookup[idx]; + let (len, s) = (u32::from(e >> 8), u32::from(e & 0xFF)); + if len == 0 { + continue; + } + let d: i32 = match s { + 0 => 0, + 16 => 32768, + 1..=15 if len + s <= LOOKUP_BITS => { + let v = ((idx as u32 >> (LOOKUP_BITS - len - s)) & ((1 << s) - 1)) as i32; + if v < (1 << (s - 1)) { v - (1 << s) + 1 } else { v } + } + _ => continue, + }; + *f = (d << 8) | (len + if (1..=15).contains(&s) { s } else { 0 }) as i32; + } + Ok(h) + } +} + +/// MSB-first bit reader over entropy-coded data: removes `FF 00` stuffing, +/// stops at markers (then feeds zero bits). +struct Bits<'a> { + data: &'a [u8], + pos: usize, + acc: u64, + n: u32, + marker: bool, + /// Zero bytes fed past the end of the data or a marker. + fake: usize, +} + +impl<'a> Bits<'a> { + fn new(data: &'a [u8]) -> Self { + Bits { data, pos: 0, acc: 0, n: 0, marker: false, fake: 0 } + } + + #[inline] + fn refill(&mut self) { + // Fast path: the next 8 bytes hold no 0xFF (no stuffing, no marker). + if !self.marker + && let Some(s) = self.data.get(self.pos..self.pos + 8) + { + let chunk = u64::from_be_bytes([s[0], s[1], s[2], s[3], s[4], s[5], s[6], s[7]]); + let inv = !chunk; + if inv.wrapping_sub(0x0101_0101_0101_0101) & !inv & 0x8080_8080_8080_8080 == 0 { + let k = (64 - self.n) / 8; // whole bytes that fit (1..=8) + if k == 0 { + return; + } + if k == 8 { + self.acc = chunk; + } else { + self.acc |= (chunk >> (64 - 8 * k)) << (64 - self.n - 8 * k); + } + self.n += 8 * k; + self.pos += k as usize; + return; + } + } + while self.n <= 56 { + let byte = if self.marker { + self.fake += 1; + 0 + } else { + match self.data.get(self.pos) { + Some(&0xFF) => match self.data.get(self.pos + 1) { + Some(0) => { + self.pos += 2; + 0xFF + } + _ => { + self.marker = true; + self.fake += 1; + 0 + } + }, + Some(&b) => { + self.pos += 1; + b + } + None => { + self.marker = true; + self.fake += 1; + 0 + } + } + }; + self.acc |= u64::from(byte) << (56 - self.n); + self.n += 8; + } + } + + #[inline] + fn peek(&self, k: u32) -> u32 { + (self.acc >> (64 - k)) as u32 + } + + #[inline] + fn consume(&mut self, k: u32) { + self.acc <<= k; + self.n -= k; + } + + /// Reads `k` (1..=16) bits. + #[inline] + fn get(&mut self, k: u32) -> u32 { + if self.n < k { + self.refill(); + } + let v = self.peek(k); + self.consume(k); + v + } + + /// Skips to the next RSTn marker and past it; resets the bit buffer. + fn restart(&mut self) -> Result<()> { + self.acc = 0; + self.n = 0; + self.marker = false; + self.fake = 0; + // Tolerate fill bytes before the marker. + while let Some(&b) = self.data.get(self.pos) { + if b == 0xFF { + match self.data.get(self.pos + 1) { + Some(0xD0..=0xD7) => { + self.pos += 2; + return Ok(()); + } + Some(0xFF) => self.pos += 1, + _ => return Err(RawError::malformed("expected a restart marker")), + } + } else { + self.pos += 1; + } + } + Err(RawError::malformed("missing restart marker")) + } + + #[inline] + fn decode(&mut self, h: &Huffman) -> Result { + if self.n < 32 { + self.refill(); + } + let e = h.lookup[self.peek(LOOKUP_BITS) as usize]; + if e != 0 { + self.consume(u32::from(e >> 8)); + return Ok(u32::from(e & 0xFF)); + } + for len in LOOKUP_BITS + 1..=16 { + let code = self.peek(len) as i32; + if code <= h.maxcode[len as usize] { + let i = h.valptr[len as usize] + (code - h.mincode[len as usize]) as usize; + let v = *h.vals.get(i).ok_or_else(|| RawError::malformed("bad Huffman code"))?; + self.consume(len); + return Ok(u32::from(v)); + } + } + Err(RawError::malformed("bad Huffman code")) + } + + /// A difference value: Huffman-coded category SSSS then SSSS extra bits (T.81 H.1.2.2). + #[inline] + fn diff(&mut self, h: &Huffman) -> Result { + if self.n < 32 { + self.refill(); + } + let f = h.fast[self.peek(LOOKUP_BITS) as usize]; + if f != 0 { + self.consume((f & 0xFF) as u32); + return Ok(f >> 8); + } + let s = self.decode(h)?; + Ok(match s { + 0 => 0, + 1..=15 => { + let v = self.get(s) as i32; + if v < (1 << (s - 1)) { v - (1 << s) + 1 } else { v } + } + 16 => 32768, + _ => return Err(RawError::malformed("bad difference category")), + }) + } +} + +fn be16(d: &[u8], at: usize) -> Result { + let s = d.get(at..at + 2).ok_or_else(|| RawError::malformed("truncated lossless JPEG"))?; + Ok(usize::from(u16::from_be_bytes([s[0], s[1]]))) +} + +/// Reads the frame header without decoding. +pub(crate) fn frame(data: &[u8]) -> Result { + parse(data).map(|p| p.frame) +} + +struct Parsed<'a> { + frame: Frame, + /// Per component of the scan (in frame order): its Huffman table. + tables: Vec, + predictor: u8, + point_transform: u8, + restart: usize, + entropy: &'a [u8], +} + +fn parse(data: &[u8]) -> Result> { + if data.get(0..2) != Some(&[0xFF, 0xD8]) { + return Err(RawError::malformed("lossless JPEG: missing SOI")); + } + let mut pos = 2; + let mut frame: Option = None; + let mut comp_ids: Vec = Vec::new(); + let mut huff: [Option; 4] = [None, None, None, None]; + let mut restart = 0usize; + loop { + // Skip fill bytes up to the marker. + while data.get(pos) == Some(&0xFF) && data.get(pos + 1) == Some(&0xFF) { + pos += 1; + } + if data.get(pos) != Some(&0xFF) { + return Err(RawError::malformed("lossless JPEG: expected a marker")); + } + let marker = *data.get(pos + 1).ok_or_else(|| RawError::malformed("truncated lossless JPEG"))?; + pos += 2; + if marker == 0xD9 { + return Err(RawError::malformed("lossless JPEG: no scan")); + } + if (0xD0..=0xD7).contains(&marker) || marker == 0x01 { + continue; + } + let len = be16(data, pos)?; + if len < 2 { + return Err(RawError::malformed("lossless JPEG: bad segment length")); + } + let seg = data.get(pos + 2..pos + len).ok_or_else(|| RawError::malformed("truncated lossless JPEG segment"))?; + match marker { + 0xC3 => { + if seg.len() < 6 { + return Err(RawError::malformed("short SOF3")); + } + let precision = seg[0]; + let height = usize::from(u16::from_be_bytes([seg[1], seg[2]])); + let width = usize::from(u16::from_be_bytes([seg[3], seg[4]])); + let nf = usize::from(seg[5]); + if !(2..=16).contains(&precision) || width == 0 || height == 0 || nf == 0 || nf > 4 || seg.len() < 6 + 3 * nf { + return Err(RawError::malformed(format!("unsupported SOF3 ({precision}-bit, {width}x{height}, {nf} components)"))); + } + comp_ids.clear(); + for c in 0..nf { + let s = &seg[6 + 3 * c..9 + 3 * c]; + if s[1] != 0x11 { + return Err(RawError::unsupported("lossless JPEG with subsampled components")); + } + comp_ids.push(s[0]); + } + frame = Some(Frame { precision, width, height, components: nf }); + } + 0xC0..=0xC2 => return Err(RawError::unsupported("lossy (DCT) JPEG compressed raw data")), + 0xC5..=0xC7 | 0xC9..=0xCB | 0xCD..=0xCF => return Err(RawError::unsupported("hierarchical or arithmetic-coded JPEG")), + 0xC4 => { + let mut p = 0; + while p < seg.len() { + let tc_th = seg[p]; + let id = usize::from(tc_th & 0x0F); + let counts: [u8; 16] = seg.get(p + 1..p + 17).and_then(|s| s.try_into().ok()).ok_or_else(|| RawError::malformed("short DHT"))?; + let total: usize = counts.iter().map(|&c| usize::from(c)).sum(); + let vals = seg.get(p + 17..p + 17 + total).ok_or_else(|| RawError::malformed("short DHT"))?; + if id < 4 { + huff[id] = Some(Huffman::new(&counts, vals)?); + } + p += 17 + total; + } + } + 0xDD => { + restart = be16(seg, 0)?; + } + 0xDA => { + let f = frame.ok_or_else(|| RawError::malformed("lossless JPEG: SOS before SOF3"))?; + let ns = usize::from(*seg.first().ok_or_else(|| RawError::malformed("short SOS"))?); + if ns != f.components || seg.len() < 1 + 2 * ns + 3 { + return Err(RawError::unsupported("lossless JPEG with non-interleaved scans")); + } + let mut tables = Vec::with_capacity(ns); + for c in 0..ns { + let id = seg[1 + 2 * c]; + if comp_ids.get(c) != Some(&id) { + return Err(RawError::unsupported("lossless JPEG scan component order")); + } + let td = usize::from(seg[2 + 2 * c] >> 4); + let h = huff.get(td).and_then(Option::as_ref).ok_or_else(|| RawError::malformed("missing Huffman table"))?; + tables.push(h.clone()); + } + let predictor = seg[1 + 2 * ns]; + let point_transform = seg[3 + 2 * ns] & 0x0F; + if !(1..=7).contains(&predictor) { + return Err(RawError::malformed(format!("lossless JPEG predictor {predictor}"))); + } + if point_transform >= f.precision { + return Err(RawError::malformed("lossless JPEG point transform")); + } + let entropy = data.get(pos + len..).unwrap_or(&[]); + return Ok(Parsed { frame: f, tables, predictor, point_transform, restart, entropy }); + } + _ => {} // APPn, COM, DQT, DNL… ignored + } + pos += len; + } +} + +/// Decodes a lossless JPEG stream. `max_samples` bounds the output size and +/// is checked before allocating. +pub(crate) fn decode(data: &[u8], max_samples: usize) -> Result<(Frame, Vec)> { + let p = parse(data)?; + let f = p.frame; + let total = f + .samples() + .filter(|&n| n <= max_samples) + .ok_or_else(|| RawError::LimitExceeded(format!("lossless JPEG frame {}x{}x{}", f.width, f.height, f.components)))?; + // An entropy-coded sample needs at least one bit: a stream far shorter than + // the frame cannot be valid (cheap decompression-bomb guard). + if p.entropy.len().saturating_mul(8) < total / 2 { + return Err(RawError::malformed("lossless JPEG data is truncated")); + } + let row_len = f.width * f.components; + if p.restart != 0 && p.restart % f.width != 0 { + return Err(RawError::unsupported("lossless JPEG restart interval that is not whole rows")); + } + let rows_per_interval = if p.restart == 0 { usize::MAX } else { p.restart / f.width }; + let mut out = vec![0u16; total]; + let mut bits = Bits::new(p.entropy); + let nc = f.components; + let initial = 1i32 << (f.precision - p.point_transform - 1); + let pred = p.predictor; + let mut interval_row = 0usize; + for y in 0..f.height { + if interval_row == rows_per_interval { + bits.restart()?; + interval_row = 0; + } + let first_line = interval_row == 0; + let (prev, cur) = if y == 0 { + let (_, cur) = out.split_at_mut(0); + (&[][..], &mut cur[..row_len]) + } else { + let (a, b) = out.split_at_mut(y * row_len); + (&a[(y - 1) * row_len..], &mut b[..row_len]) + }; + if pred == 1 || first_line { + // Predictor 1 (or the first line of an interval): Ra after the first column. + for c in 0..nc { + let px = if first_line { initial } else { i32::from(prev[c]) }; + cur[c] = (px + bits.diff(&p.tables[c])?) as u16; + } + if nc == 1 { + let t = &p.tables[0]; + for i in 1..row_len { + cur[i] = (i32::from(cur[i - 1]) + bits.diff(t)?) as u16; + } + } else { + for i in (nc..row_len).step_by(nc) { + for c in 0..nc { + cur[i + c] = (i32::from(cur[i + c - nc]) + bits.diff(&p.tables[c])?) as u16; + } + } + } + interval_row += 1; + if bits.fake > 64 { + return Err(RawError::malformed("lossless JPEG data is truncated")); + } + continue; + } + for x in 0..f.width { + for c in 0..nc { + let i = x * nc + c; + let px = if x == 0 { + if first_line { initial } else { i32::from(prev[i]) } + } else { + let ra = i32::from(cur[i - nc]); + if first_line { + ra + } else { + let rb = i32::from(prev[i]); + let rc = i32::from(prev[i - nc]); + match pred { + 1 => ra, + 2 => rb, + 3 => rc, + 4 => ra + rb - rc, + 5 => ra + ((rb - rc) >> 1), + 6 => rb + ((ra - rc) >> 1), + _ => (ra + rb) >> 1, + } + } + }; + let d = bits.diff(&p.tables[c])?; + cur[i] = (px + d) as u16; + } + } + interval_row += 1; + if bits.fake > 64 { + return Err(RawError::malformed("lossless JPEG data is truncated")); + } + } + if p.point_transform > 0 { + for v in &mut out { + *v <<= p.point_transform; + } + } + Ok((f, out)) +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::testgen::lj92_encode; + + /// Deterministic pseudo-random samples with some structure (smooth + noise). + fn samples(n: usize, precision: u8, seed: u64) -> Vec { + let mut s = seed; + let max = (1u32 << precision) - 1; + (0..n) + .map(|i| { + s = s.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407); + let noise = (s >> 33) as u32 % 64; + let base = (i as u32 * 37) % (max + 1); + ((base + noise) % (max + 1)) as u16 + }) + .collect() + } + + #[test] + fn roundtrip_all_predictors_and_precisions() { + for precision in [2u8, 8, 12, 14, 16] { + for predictor in 1..=7u8 { + for comps in 1..=4usize { + let (w, h) = (13, 7); + let v = samples(w * h * comps, precision, u64::from(predictor) * 31 + comps as u64); + let j = lj92_encode(&v, w, h, comps, precision, predictor, None); + let (f, out) = decode(&j, usize::MAX).unwrap(); + assert_eq!(f, Frame { precision, width: w, height: h, components: comps }); + assert_eq!(out, v, "P={precision} pred={predictor} comps={comps}"); + } + } + } + } + + #[test] + fn roundtrip_with_restarts_and_extreme_values() { + let (w, h) = (16, 9); + // Alternating 0 / 65535 forces the 16-bit and 32768 difference categories. + let v: Vec = (0..w * h * 2).map(|i| if (i / 3) % 2 == 0 { 0 } else { 65535 }).collect(); + for rows in [1, 2, 4] { + let j = lj92_encode(&v, w, h, 2, 16, 1, Some(rows)); + let (_, out) = decode(&j, usize::MAX).unwrap(); + assert_eq!(out, v, "restart every {rows} rows"); + } + } + + #[test] + fn hostile_streams_error_without_panicking() { + let v = samples(32 * 8, 14, 7); + let j = lj92_encode(&v, 32, 8, 1, 14, 1, Some(2)); + for n in 0..j.len() { + let _ = decode(&j[..n], usize::MAX); + } + let mut s = 0x1234_5678u64; + for _ in 0..2000 { + let mut b = j.clone(); + for _ in 0..4 { + s = s.wrapping_mul(6364136223846793005).wrapping_add(1); + let at = (s >> 33) as usize % b.len(); + b[at] = (s >> 20) as u8; + } + let _ = decode(&b, 1 << 20); + } + // A frame larger than allowed is refused before allocating. + assert!(matches!(decode(&j, 10), Err(RawError::LimitExceeded(_)))); + // A huge frame with a tiny entropy segment is refused as truncated. + let mut big = j.clone(); + big[7..9].copy_from_slice(&60000u16.to_be_bytes()); // height + big[9..11].copy_from_slice(&60000u16.to_be_bytes()); // width + assert!(decode(&big, usize::MAX).is_err()); + // Lossy JPEG markers are reported as unsupported. + let lossy = [0xFF, 0xD8, 0xFF, 0xC0, 0x00, 0x0B, 8, 0, 1, 0, 1, 1, 1, 0x11, 0]; + assert!(matches!(decode(&lossy, usize::MAX), Err(RawError::Unsupported(_)))); + } +} diff --git a/crates/raw/src/opcodes.rs b/crates/raw/src/opcodes.rs new file mode 100644 index 00000000..9dc231ca --- /dev/null +++ b/crates/raw/src/opcodes.rs @@ -0,0 +1,219 @@ +//! DNG opcode lists (DNG 1.7, chapter 7 "Opcode List Processing"). +//! +//! An opcode list is big-endian whatever the file's byte order: a count, then +//! per opcode its ID, DNG version, flags (bit 0 = optional), parameter size +//! and parameters. Only GainMap (ID 9), the lens-shading correction phones +//! write in OpcodeList2, is applied; other opcodes are reported. + +/// One GainMap opcode: a grid of gains, bilinearly interpolated, applied to +/// the pixels of an area (in active-area coordinates) at a row / column pitch. +#[derive(Debug, Clone, PartialEq)] +pub struct GainMap { + pub top: usize, + pub left: usize, + pub bottom: usize, + pub right: usize, + pub plane: usize, + pub planes: usize, + pub row_pitch: usize, + pub col_pitch: usize, + pub points_v: usize, + pub points_h: usize, + /// Spacing and origin of the map points, in image-relative coordinates (0..1). + pub spacing_v: f64, + pub spacing_h: f64, + pub origin_v: f64, + pub origin_h: f64, + pub map_planes: usize, + /// `points_v × points_h × map_planes` gains. + pub gains: Vec, +} + +impl GainMap { + /// Gains along the map's horizontal points for image-relative row `rv` + /// and map plane `mp` (vertical interpolation done). + pub fn row(&self, rv: f64, mp: usize) -> Vec { + let mp = mp.min(self.map_planes.saturating_sub(1)); + let fv = ((rv - self.origin_v) / self.spacing_v).clamp(0.0, (self.points_v - 1) as f64); + let i0 = (fv.floor() as usize).min(self.points_v - 1); + let i1 = (i0 + 1).min(self.points_v - 1); + let t = (fv - i0 as f64) as f32; + (0..self.points_h) + .map(|j| { + let a = self.gains.get((i0 * self.points_h + j) * self.map_planes + mp).copied().unwrap_or(1.0); + let b = self.gains.get((i1 * self.points_h + j) * self.map_planes + mp).copied().unwrap_or(1.0); + a + (b - a) * t + }) + .collect() + } + + /// Interpolates a [`GainMap::row`] at image-relative column `rh`. + #[inline] + pub fn at(&self, row: &[f32], rh: f64) -> f32 { + let fh = ((rh - self.origin_h) / self.spacing_h).clamp(0.0, (self.points_h - 1) as f64); + let j0 = (fh.floor() as usize).min(self.points_h - 1); + let j1 = (j0 + 1).min(self.points_h - 1); + let t = (fh - j0 as f64) as f32; + let (a, b) = (row.get(j0).copied().unwrap_or(1.0), row.get(j1).copied().unwrap_or(1.0)); + a + (b - a) * t + } + + /// `true` if the map covers active-area pixel (`ax`, `ay`). + #[inline] + pub fn covers(&self, ax: usize, ay: usize) -> bool { + ay >= self.top + && ay < self.bottom + && ax >= self.left + && ax < self.right + && (ay - self.top).is_multiple_of(self.row_pitch) + && (ax - self.left).is_multiple_of(self.col_pitch) + } +} + +/// Opcode names (DNG 1.7 table) for warnings. +fn name(id: u32) -> String { + match id { + 1 => "WarpRectilinear".into(), + 2 => "WarpFisheye".into(), + 3 => "FixVignetteRadial".into(), + 4 => "FixBadPixelsConstant".into(), + 5 => "FixBadPixelsList".into(), + 6 => "TrimBounds".into(), + 7 => "MapTable".into(), + 8 => "MapPolynomial".into(), + 9 => "GainMap".into(), + 10 => "DeltaPerRow".into(), + 11 => "DeltaPerColumn".into(), + 12 => "ScalePerRow".into(), + 13 => "ScalePerColumn".into(), + 14 => "WarpRectilinear2".into(), + id => format!("opcode {id}"), + } +} + +fn be32(b: &[u8], at: usize) -> Option { + let s: [u8; 4] = b.get(at..at.checked_add(4)?)?.try_into().ok()?; + Some(u32::from_be_bytes(s)) +} + +fn be64f(b: &[u8], at: usize) -> Option { + let s: [u8; 8] = b.get(at..at.checked_add(8)?)?.try_into().ok()?; + Some(f64::from_be_bytes(s)).filter(|v| v.is_finite()) +} + +fn gain_map(p: &[u8]) -> Option { + let u = |i: usize| be32(p, 4 * i).map(|v| v as usize); + let (top, left, bottom, right) = (u(0)?, u(1)?, u(2)?, u(3)?); + let (plane, planes, row_pitch, col_pitch) = (u(4)?, u(5)?, u(6)?, u(7)?); + let (points_v, points_h) = (u(8)?, u(9)?); + let (spacing_v, spacing_h, origin_v, origin_h) = (be64f(p, 40)?, be64f(p, 48)?, be64f(p, 56)?, be64f(p, 64)?); + let map_planes = u(18)?; + let n = points_v.checked_mul(points_h)?.checked_mul(map_planes)?; + if n == 0 || n > 1 << 20 || row_pitch == 0 || col_pitch == 0 || spacing_v <= 0.0 || spacing_h <= 0.0 || bottom <= top || right <= left { + return None; + } + let gains = (0..n).map(|i| be32(p, 76 + 4 * i).map(f32::from_bits).filter(|g| g.is_finite())).collect::>>()?; + Some(GainMap { + top, + left, + bottom, + right, + plane, + planes, + row_pitch, + col_pitch, + points_v, + points_h, + spacing_v, + spacing_h, + origin_v, + origin_h, + map_planes, + gains, + }) +} + +/// Parses an opcode list: the gain maps, and the names of opcodes not applied. +pub(crate) fn parse(list: &[u8]) -> (Vec, Vec) { + let mut maps = Vec::new(); + let mut skipped = Vec::new(); + let Some(count) = be32(list, 0) else { return (maps, skipped) }; + let mut at = 4usize; + for _ in 0..count.min(4096) { + let (Some(id), Some(size)) = (be32(list, at), be32(list, at.saturating_add(12))) else { break }; + let start = at.saturating_add(16); + let Some(params) = list.get(start..start.saturating_add(size as usize)) else { break }; + match (id, gain_map(params)) { + (9, Some(m)) => maps.push(m), + _ => { + let n = name(id); + if !skipped.contains(&n) { + skipped.push(n); + } + } + } + at = at.saturating_add(16).saturating_add(size as usize); + } + (maps, skipped) +} + +#[cfg(test)] +mod tests { + use super::*; + + fn list(ops: &[(u32, Vec)]) -> Vec { + let mut b = (ops.len() as u32).to_be_bytes().to_vec(); + for (id, p) in ops { + b.extend_from_slice(&id.to_be_bytes()); + b.extend_from_slice(&0x0103_0000u32.to_be_bytes()); + b.extend_from_slice(&1u32.to_be_bytes()); + b.extend_from_slice(&(p.len() as u32).to_be_bytes()); + b.extend_from_slice(p); + } + b + } + + pub(crate) fn gain_map_params(area: [u32; 4], pitch: u32, gains: &[[f32; 2]; 2]) -> Vec { + let mut p = Vec::new(); + for v in [area[0], area[1], area[2], area[3], 0, 1, pitch, pitch, 2, 2] { + p.extend_from_slice(&v.to_be_bytes()); + } + for v in [1.0f64, 1.0, 0.0, 0.0] { + p.extend_from_slice(&v.to_be_bytes()); + } + p.extend_from_slice(&1u32.to_be_bytes()); + for row in gains { + for g in row { + p.extend_from_slice(&g.to_bits().to_be_bytes()); + } + } + p + } + + #[test] + fn parses_gain_maps_and_names_others() { + let b = list(&[(9, gain_map_params([0, 0, 10, 10], 1, &[[1.0, 2.0], [3.0, 4.0]])), (1, vec![0; 8])]); + let (maps, skipped) = parse(&b); + assert_eq!(maps.len(), 1); + assert_eq!(skipped, vec!["WarpRectilinear".to_string()]); + let m = &maps[0]; + // Bilinear: centre of the 2×2 grid. + let row = m.row(0.5, 0); + assert!((m.at(&row, 0.5) - 2.5).abs() < 1e-6); + assert!((m.at(&m.row(0.0, 0), 0.0) - 1.0).abs() < 1e-6); + assert!((m.at(&m.row(1.0, 0), 1.0) - 4.0).abs() < 1e-6); + // Outside the grid clamps. + assert!((m.at(&m.row(-3.0, 0), 9.0) - 2.0).abs() < 1e-6); + } + + #[test] + fn hostile_lists_do_not_panic() { + let b = list(&[(9, gain_map_params([0, 0, 10, 10], 2, &[[1.0, 2.0], [3.0, 4.0]]))]); + for n in 0..b.len() { + let _ = parse(&b[..n]); + } + let mut huge = b.clone(); + huge[4 + 16 + 32..4 + 16 + 40].copy_from_slice(&[0xFF; 8]); // points + assert!(parse(&huge).0.is_empty()); + } +} diff --git a/crates/raw/src/par.rs b/crates/raw/src/par.rs new file mode 100644 index 00000000..a8e2c9f9 --- /dev/null +++ b/crates/raw/src/par.rs @@ -0,0 +1,29 @@ +//! Data-parallel helpers: rayon on native targets, sequential on wasm. + +/// Calls `f(index, chunk)` for each `chunk`-sized piece of `buf`. +pub(crate) fn chunks_mut(buf: &mut [T], chunk: usize, f: impl Fn(usize, &mut [T]) + Sync + Send) { + let chunk = chunk.max(1); + #[cfg(not(target_arch = "wasm32"))] + { + use rayon::prelude::*; + buf.par_chunks_mut(chunk).enumerate().for_each(|(i, c)| f(i, c)); + } + #[cfg(target_arch = "wasm32")] + buf.chunks_mut(chunk).enumerate().for_each(|(i, c)| f(i, c)); +} + +/// Maps `f` over `0..n`, keeping order. +pub(crate) fn map(n: usize, f: impl Fn(usize) -> T + Sync + Send) -> Vec { + #[cfg(not(target_arch = "wasm32"))] + { + use rayon::prelude::*; + (0..n).into_par_iter().map(f).collect() + } + #[cfg(target_arch = "wasm32")] + (0..n).map(f).collect() +} + +/// Rows per parallel band for an image `width` samples wide (about 64 K samples a band). +pub(crate) fn band_rows(width: usize) -> usize { + (65536 / width.max(1)).clamp(1, 256) +} diff --git a/crates/raw/src/preview.rs b/crates/raw/src/preview.rs new file mode 100644 index 00000000..adeb6865 --- /dev/null +++ b/crates/raw/src/preview.rs @@ -0,0 +1,81 @@ +//! Finding the camera-rendered JPEG preview embedded in a raw file. + +use crate::tiff::{Tiff, tag}; + +/// An embedded baseline-JPEG preview. +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub struct Preview<'a> { + /// The complete JPEG stream. + pub jpeg: &'a [u8], + pub width: u32, + pub height: u32, +} + +/// Dimensions of a JPEG whose frame is DCT-coded (baseline / extended / +/// progressive), or `None` (lossless JPEG raw data, or not a JPEG). +fn dct_jpeg_size(b: &[u8]) -> Option<(u32, u32)> { + if b.get(0..3) != Some(&[0xFF, 0xD8, 0xFF]) { + return None; + } + let mut pos = 2usize; + for _ in 0..256 { + if *b.get(pos)? != 0xFF { + return None; + } + let m = *b.get(pos + 1)?; + if m == 0xFF { + pos += 1; + continue; + } + let len = usize::from(u16::from_be_bytes([*b.get(pos + 2)?, *b.get(pos + 3)?])); + if matches!(m, 0xC0..=0xC2) { + let h = u16::from_be_bytes([*b.get(pos + 5)?, *b.get(pos + 6)?]); + let w = u16::from_be_bytes([*b.get(pos + 7)?, *b.get(pos + 8)?]); + return (w > 0 && h > 0).then_some((u32::from(w), u32::from(h))); + } + if matches!(m, 0xC3 | 0xC5..=0xC7 | 0xC9..=0xCB | 0xCD..=0xCF | 0xDA | 0xD9) || len < 2 { + return None; + } + pos = pos.checked_add(2 + len)?; + } + None +} + +fn candidate(data: &[u8], offset: usize, len: usize) -> Option> { + let jpeg = data.get(offset..offset.checked_add(len)?)?; + let (width, height) = dct_jpeg_size(jpeg)?; + Some(Preview { jpeg, width, height }) +} + +/// The largest embedded JPEG preview of a raw file (TIFF-based raws and RAF). +pub fn embedded_preview(bytes: &[u8]) -> Option> { + let mut found: Vec = Vec::new(); + if bytes.starts_with(b"FUJIFILMCCD-RAW") { + // RAF header: big-endian JPEG offset and length at bytes 84 and 88. + let be = |o: usize| bytes.get(o..o + 4).map(|s| u32::from_be_bytes([s[0], s[1], s[2], s[3]]) as usize); + if let (Some(off), Some(len)) = (be(84), be(88)) { + found.extend(candidate(bytes, off, len)); + } + } else if let Some(t) = Tiff::new(bytes) { + for ifd in t.all_ifds() { + if let (Some(off), Some(len)) = (t.tag_uint(&ifd, tag::JPEG_INTERCHANGE_FORMAT), t.tag_uint(&ifd, tag::JPEG_INTERCHANGE_FORMAT_LENGTH)) { + found.extend(candidate(bytes, off as usize, len as usize)); + } + if matches!(t.tag_uint(&ifd, tag::COMPRESSION), Some(6 | 7)) { + let offs = t.tag_uints(&ifd, tag::STRIP_OFFSETS); + let lens = t.tag_uints(&ifd, tag::STRIP_BYTE_COUNTS); + if let ([off], [len]) = (offs.as_slice(), lens.as_slice()) { + found.extend(candidate(bytes, *off as usize, *len as usize)); + } + } + // Panasonic RW2 JpgFromRaw (0x002E) holds a complete JPEG. + if let Some(e) = ifd.get(0x002E).filter(|e| e.typ == 7) + && let Some(raw) = t.raw(e) + && let Some((width, height)) = dct_jpeg_size(raw) + { + found.push(Preview { jpeg: raw, width, height }); + } + } + } + found.into_iter().max_by_key(|p| u64::from(p.width) * u64::from(p.height)) +} diff --git a/crates/raw/src/sensor.rs b/crates/raw/src/sensor.rs new file mode 100644 index 00000000..59ec38f5 --- /dev/null +++ b/crates/raw/src/sensor.rs @@ -0,0 +1,366 @@ +//! Undeveloped sensor data and the shared strip/tile reader. + +use crate::color::ColorInfo; +use crate::error::{RawError, Result}; +use crate::tiff::{Ifd, Tiff, tag}; +use crate::{Limits, RawFormat, ljpeg, par}; + +/// A rectangle in sensor-data pixels. +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub struct Rect { + pub x: usize, + pub y: usize, + pub width: usize, + pub height: usize, +} + +impl Rect { + pub fn new(x: usize, y: usize, width: usize, height: usize) -> Self { + Rect { x, y, width, height } + } + + /// The overlap of two rectangles (may be empty). + pub fn intersect(&self, o: &Rect) -> Rect { + let x0 = self.x.max(o.x); + let y0 = self.y.max(o.y); + let x1 = (self.x.saturating_add(self.width)).min(o.x.saturating_add(o.width)); + let y1 = (self.y.saturating_add(self.height)).min(o.y.saturating_add(o.height)); + Rect::new(x0, y0, x1.saturating_sub(x0), y1.saturating_sub(y0)) + } + + pub fn is_empty(&self) -> bool { + self.width == 0 || self.height == 0 + } +} + +/// Colour of one CFA cell: 0 = red, 1 = green, 2 = blue. +#[derive(Debug, Clone, PartialEq, Eq)] +pub struct Cfa { + pub width: usize, + pub height: usize, + /// Row-major colours, `width * height` entries. + pub colors: Vec, + /// Data coordinates of the pattern's top-left cell. + pub origin_x: usize, + pub origin_y: usize, +} + +impl Cfa { + /// A 2×2 pattern from a string such as `"RGGB"`, anchored at (0, 0). + pub fn bayer(pattern: &str) -> Option { + let colors: Vec = pattern + .chars() + .map(|c| match c { + 'R' => Some(0), + 'G' => Some(1), + 'B' => Some(2), + _ => None, + }) + .collect::>()?; + (colors.len() == 4).then_some(Cfa { width: 2, height: 2, colors, origin_x: 0, origin_y: 0 }) + } + + /// Colour at data coordinates. + #[inline] + pub fn color(&self, x: usize, y: usize) -> u8 { + let cx = (x + self.width - self.origin_x % self.width) % self.width; + let cy = (y + self.height - self.origin_y % self.height) % self.height; + self.colors.get(cy * self.width + cx).copied().unwrap_or(1) + } + + /// The 2×2 colours seen from data coordinates (`x0`, `y0`), row-major. + pub fn phase(&self, x0: usize, y0: usize) -> [u8; 4] { + [self.color(x0, y0), self.color(x0 + 1, y0), self.color(x0, y0 + 1), self.color(x0 + 1, y0 + 1)] + } + + /// `true` for a Bayer pattern: a 2×2 period holding one red, one blue and two + /// diagonal greens. + pub fn is_bayer(&self) -> bool { + if self.width == 0 || self.height == 0 || self.width > 16 || self.height > 16 { + return false; + } + let p = self.phase(self.origin_x, self.origin_y); + let periodic = (0..self.height).all(|y| (0..self.width).all(|x| self.color(self.origin_x + x, self.origin_y + y) == p[(y % 2) * 2 + x % 2])); + let greens_on_diagonal = + (p[0] == 1 && p[3] == 1 && p[1] != 1 && p[2] != 1 && p[1] + p[2] == 2) || (p[1] == 1 && p[2] == 1 && p[0] != 1 && p[3] != 1 && p[0] + p[3] == 2); + periodic && greens_on_diagonal + } +} + +/// Black levels (DNG BlackLevel / BlackLevelRepeatDim / BlackLevelDeltaH/V). +#[derive(Debug, Clone, PartialEq)] +pub struct BlackLevels { + pub rows: usize, + pub cols: usize, + /// `rows * cols * samples` values, sample fastest. + pub values: Vec, + /// Per column of the active area. + pub delta_h: Vec, + /// Per row of the active area. + pub delta_v: Vec, +} + +impl BlackLevels { + pub fn uniform(v: f32) -> Self { + BlackLevels { rows: 1, cols: 1, values: vec![v], delta_h: Vec::new(), delta_v: Vec::new() } + } + + /// Black level at active-area coordinates for sample `s` of `spp`. + #[inline] + pub fn at(&self, ax: usize, ay: usize, s: usize, spp: usize) -> f32 { + let r = ay % self.rows.max(1); + let c = ax % self.cols.max(1); + let base = self.values.get((r * self.cols.max(1) + c) * spp + s).or_else(|| self.values.first()).copied().unwrap_or(0.0); + base + self.delta_h.get(ax).copied().unwrap_or(0.0) + self.delta_v.get(ay).copied().unwrap_or(0.0) + } + + /// `true` when there is a single black value for everything. + pub fn is_uniform(&self) -> bool { + self.delta_h.iter().all(|v| *v == 0.0) && self.delta_v.iter().all(|v| *v == 0.0) && self.values.windows(2).all(|w| w[0] == w[1]) + } +} + +/// Undeveloped sensor data with everything needed to develop it. +#[derive(Debug, Clone)] +pub struct Sensor { + pub format: RawFormat, + pub make: Option, + pub model: Option, + /// Data width and height in pixels. + pub width: usize, + pub height: usize, + /// Samples per pixel: 1 for CFA data, 3 for LinearRaw. + pub samples: usize, + /// Row-major interleaved samples, `width * height * samples`. + pub data: Vec, + /// The CFA pattern (CFA data only). + pub cfa: Option, + /// Maps stored values to linear values (DNG LinearizationTable). + pub linearization: Option>, + pub black: BlackLevels, + /// White (saturation) level per sample. + pub white: [f32; 3], + /// Sensor area holding image data (DNG ActiveArea); black-level deltas are indexed from it. + pub active: Rect, + /// The area to output (ActiveArea ∩ DefaultCrop). + pub crop: Rect, + pub color: ColorInfo, + /// White-balance multipliers from the camera's own metadata (non-DNG), R G B. + pub camera_wb: Option<[f64; 3]>, + /// TIFF orientation (1–8). + pub orientation: u16, + /// Exposure compensation baked into the default rendering (DNG BaselineExposure), in EV. + pub baseline_exposure: f64, + /// DNG OpcodeList2 gain maps (lens shading), applied to the normalized data. + pub gain_maps: Vec, + pub warnings: Vec, +} + +/// Raw samples read from one IFD. +pub(crate) struct Plane { + pub width: usize, + pub height: usize, + pub samples: usize, + pub bits: u32, + pub data: Vec, +} + +struct Segment { + x: usize, + y: usize, + /// Encoded width and height (tile size, or image width × rows per strip). + w: usize, + h: usize, + offset: usize, + len: usize, +} + +/// How a lossless-JPEG tile's samples map onto the tile. +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub(crate) enum JpegLayout { + /// The decoded sample sequence fills the tile row by row (DNG). + Flat, + /// A 4-component frame of half the tile's width and height: each JPEG + /// pixel holds a 2×2 block, row-major (observed in Sony lossless ARW). + Quads, +} + +/// Reads the (uncompressed or lossless-JPEG) image of a raw IFD. +pub(crate) fn read_plane(t: &Tiff, ifd: &Ifd, limits: &Limits, layout: JpegLayout) -> Result { + let width = t.tag_uint(ifd, tag::IMAGE_WIDTH).ok_or_else(|| RawError::malformed("raw image has no width"))? as usize; + let height = t.tag_uint(ifd, tag::IMAGE_LENGTH).ok_or_else(|| RawError::malformed("raw image has no height"))? as usize; + let samples = t.tag_uint(ifd, tag::SAMPLES_PER_PIXEL).unwrap_or(1) as usize; + let bits = t.tag_uint(ifd, tag::BITS_PER_SAMPLE).unwrap_or(1); + let compression = t.tag_uint(ifd, tag::COMPRESSION).unwrap_or(1); + if !(1..=4).contains(&samples) { + return Err(RawError::unsupported(format!("{samples} samples per pixel"))); + } + if t.tag_uint(ifd, tag::SAMPLE_FORMAT) == Some(3) { + return Err(RawError::unsupported("floating-point raw data")); + } + if !(1..=16).contains(&bits) { + return Err(RawError::unsupported(format!("{bits}-bit raw samples"))); + } + if samples > 1 && t.tag_uint(ifd, tag::PLANAR_CONFIGURATION) == Some(2) { + return Err(RawError::unsupported("planar raw data")); + } + match compression { + 1 | 7 => {} + 8 => return Err(RawError::unsupported("Deflate-compressed (floating-point) DNG")), + 34892 => return Err(RawError::unsupported("lossy-compressed DNG")), + 52546 => return Err(RawError::unsupported("JPEG XL-compressed DNG")), + 34713 => return Err(RawError::unsupported("Nikon compressed NEF")), + 32767 => return Err(RawError::unsupported("Sony compressed ARW")), + 65535 => return Err(RawError::unsupported("Pentax compressed PEF")), + c => return Err(RawError::unsupported(format!("raw compression scheme {c}"))), + } + limits.check(width as u64, height as u64, samples as u64 * 2)?; + let segments = segments(t, ifd, width, height)?; + let le = t.le; + let decoded: Vec>> = par::map(segments.len(), |i| { + let s = &segments[i]; + let rows = s.h.min(height - s.y); + let src = t.bytes(s.offset, s.len).ok_or_else(|| RawError::malformed("raw data lies outside the file"))?; + match compression { + 1 => unpack(src, s.w, rows, samples, bits, le), + _ => { + let need = s.w * rows * samples; + let max = s.w.saturating_mul(s.h).saturating_mul(samples).saturating_mul(4); + let (f, v) = ljpeg::decode(src, max)?; + if v.len() < need { + return Err(RawError::malformed("lossless JPEG tile is smaller than its tile")); + } + if layout == JpegLayout::Quads && samples == 1 && f.components == 4 && f.width * 2 == s.w && f.height * 2 == s.h { + let mut tile = vec![0u16; s.w * s.h]; + for (i, q) in v.as_chunks::<4>().0.iter().enumerate() { + let at = 2 * (i / f.width) * s.w + 2 * (i % f.width); + tile[at] = q[0]; + tile[at + 1] = q[1]; + tile[at + s.w] = q[2]; + tile[at + s.w + 1] = q[3]; + } + return Ok(tile); + } + Ok(v) + } + } + }); + let mut data = vec![0u16; width * height * samples]; + for (s, d) in segments.iter().zip(decoded) { + let d = d?; + let rows = s.h.min(height - s.y); + let cols = s.w.min(width - s.x); + for r in 0..rows { + let src = &d[r * s.w * samples..(r * s.w + cols) * samples]; + let at = ((s.y + r) * width + s.x) * samples; + data[at..at + cols * samples].copy_from_slice(src); + } + } + Ok(Plane { width, height, samples, bits, data }) +} + +fn segments(t: &Tiff, ifd: &Ifd, width: usize, height: usize) -> Result> { + let mut out = Vec::new(); + if let (Some(tw), Some(th)) = (t.tag_uint(ifd, tag::TILE_WIDTH), t.tag_uint(ifd, tag::TILE_LENGTH)) { + let (tw, th) = (tw as usize, th as usize); + if tw == 0 || th == 0 || tw > 1 << 16 || th > 1 << 16 { + return Err(RawError::malformed("bad tile size")); + } + let across = width.div_ceil(tw); + let down = height.div_ceil(th); + let offsets = t.tag_uints(ifd, tag::TILE_OFFSETS); + let counts = t.tag_uints(ifd, tag::TILE_BYTE_COUNTS); + let n = across.checked_mul(down).ok_or_else(|| RawError::malformed("too many tiles"))?; + if offsets.len() < n || counts.len() < n { + return Err(RawError::malformed("missing tile offsets")); + } + for i in 0..n { + out.push(Segment { x: (i % across) * tw, y: (i / across) * th, w: tw, h: th, offset: offsets[i] as usize, len: counts[i] as usize }); + } + } else { + let rps = match t.tag_uint(ifd, tag::ROWS_PER_STRIP) { + Some(r) if r > 0 => (r as usize).min(height), + _ => height, + }; + let n = height.div_ceil(rps.max(1)); + let offsets = t.tag_uints(ifd, tag::STRIP_OFFSETS); + let mut counts = t.tag_uints(ifd, tag::STRIP_BYTE_COUNTS); + if offsets.len() < n { + return Err(RawError::malformed("missing strip offsets")); + } + if counts.len() < n { + if n == 1 { + counts = vec![(t.data.len().saturating_sub(offsets[0] as usize)) as u32]; + } else { + return Err(RawError::malformed("missing strip byte counts")); + } + } + for i in 0..n { + out.push(Segment { x: 0, y: i * rps, w: width, h: rps, offset: offsets[i] as usize, len: counts[i] as usize }); + } + } + Ok(out) +} + +/// Unpacks uncompressed samples: 8-bit, 16-bit in the file's byte order, or +/// MSB-first bit-packed rows padded to a byte (TIFF 6.0). 9–15-bit data +/// stored in 16-bit containers (common in camera raws) is recognised by its +/// byte count. +fn unpack(src: &[u8], w: usize, rows: usize, samples: usize, bits: u32, le: bool) -> Result> { + let n = w * rows * samples; + let truncated = || RawError::malformed("uncompressed raw data is truncated"); + if bits == 16 || (bits > 8 && src.len() >= n * 2) { + let b = src.get(..n * 2).ok_or_else(truncated)?; + return Ok(b.as_chunks::<2>().0.iter().map(|c| if le { u16::from_le_bytes([c[0], c[1]]) } else { u16::from_be_bytes([c[0], c[1]]) }).collect()); + } + if bits == 8 { + return Ok(src.get(..n).ok_or_else(truncated)?.iter().map(|&v| u16::from(v)).collect()); + } + let row_bytes = (w * samples * bits as usize).div_ceil(8); + if src.len() < row_bytes * rows { + return Err(truncated()); + } + let mut out = Vec::with_capacity(n); + for r in 0..rows { + let row = &src[r * row_bytes..(r + 1) * row_bytes]; + let mut acc: u32 = 0; + let mut have = 0u32; + let mut it = row.iter(); + for _ in 0..w * samples { + while have < bits { + acc = (acc << 8) | u32::from(*it.next().unwrap_or(&0)); + have += 8; + } + have -= bits; + out.push(((acc >> have) & ((1 << bits) - 1)) as u16); + } + } + Ok(out) +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn bayer_patterns() { + for p in ["RGGB", "BGGR", "GRBG", "GBRG"] { + assert!(Cfa::bayer(p).unwrap().is_bayer(), "{p}"); + } + for p in ["RGBG", "RRGB", "GGRB"] { + assert!(!Cfa::bayer(p).unwrap().is_bayer(), "{p}"); + } + let mut c = Cfa::bayer("RGGB").unwrap(); + c.origin_x = 1; + assert_eq!(c.phase(0, 0), [1, 0, 2, 1]); // GRBG + } + + #[test] + fn unpack_12_bit() { + // Two 12-bit samples 0xABC, 0x123 packed MSB first. + assert_eq!(unpack(&[0xAB, 0xC1, 0x23], 2, 1, 1, 12, true).unwrap(), vec![0xABC, 0x123]); + // 12-bit data in 16-bit containers. + assert_eq!(unpack(&[0xBC, 0x0A, 0x23, 0x01], 2, 1, 1, 12, true).unwrap(), vec![0xABC, 0x123]); + assert!(unpack(&[0xAB], 2, 1, 1, 12, true).is_err()); + } +} diff --git a/crates/raw/src/testgen.rs b/crates/raw/src/testgen.rs new file mode 100644 index 00000000..4e229d10 --- /dev/null +++ b/crates/raw/src/testgen.rs @@ -0,0 +1,782 @@ +//! Synthetic raw-file generators for tests and fuzz seeds: a lossless-JPEG +//! (T.81 process 14) encoder, a small TIFF writer, and DNG / CR2 / TIFF-EP +//! builders. The files are written from the public specifications only and +//! are meant to exercise the decoder, not to be camera-accurate. + +#![allow(clippy::unwrap_used, clippy::expect_used, clippy::panic, clippy::indexing_slicing)] + +// ---------------------------------------------------------------- lossless JPEG + +/// Huffman code lengths assigned to categories in decreasing frequency +/// order. Kraft sum < 1 and no code is all ones. +const LENGTHS: [u8; 17] = [2, 3, 3, 3, 4, 4, 5, 5, 6, 6, 7, 8, 9, 10, 11, 12, 13]; + +fn category(d: i32) -> u32 { + if d == 0 { + 0 + } else if d == 32768 { + 16 + } else { + 32 - d.unsigned_abs().leading_zeros() + } +} + +struct BitWriter { + out: Vec, + acc: u32, + n: u32, +} + +impl BitWriter { + fn put(&mut self, v: u32, bits: u32) { + for i in (0..bits).rev() { + self.acc = (self.acc << 1) | ((v >> i) & 1); + self.n += 1; + if self.n == 8 { + let b = self.acc as u8; + self.out.push(b); + if b == 0xFF { + self.out.push(0); + } + self.acc = 0; + self.n = 0; + } + } + } + fn flush(&mut self) { + if self.n > 0 { + let pad = 8 - self.n; + self.put((1 << pad) - 1, pad); + } + } +} + +/// Encodes `samples` (raster order, `components` interleaved, values below +/// `2^precision`) as a lossless JPEG with the given predictor (1–7) and an +/// optional restart interval in rows. +pub fn lj92_encode(samples: &[u16], width: usize, height: usize, components: usize, precision: u8, predictor: u8, restart_rows: Option) -> Vec { + assert_eq!(samples.len(), width * height * components); + let row = width * components; + let interval = restart_rows.unwrap_or(usize::MAX); + // Differences per sample, following T.81 H.1.2.1 prediction rules. + let mut diffs = Vec::with_capacity(samples.len()); + let initial = 1i32 << (precision - 1); + for y in 0..height { + let first_line = y % interval == 0; + for x in 0..width { + for c in 0..components { + let i = y * row + x * components + c; + let px = if x == 0 { + if first_line { initial } else { i32::from(samples[i - row]) } + } else if first_line { + i32::from(samples[i - components]) + } else { + let ra = i32::from(samples[i - components]); + let rb = i32::from(samples[i - row]); + let rc = i32::from(samples[i - row - components]); + match predictor { + 1 => ra, + 2 => rb, + 3 => rc, + 4 => ra + rb - rc, + 5 => ra + ((rb - rc) >> 1), + 6 => rb + ((ra - rc) >> 1), + _ => (ra + rb) >> 1, + } + }; + let mut d = (i32::from(samples[i]) - px) & 0xFFFF; + if d > 32768 { + d -= 65536; + } + diffs.push(d); + } + } + } + // Table: categories ordered by frequency get increasing code lengths. + let mut freq = [0usize; 17]; + for &d in &diffs { + freq[category(d) as usize] += 1; + } + let mut order: Vec = (0..17).collect(); + order.sort_by_key(|&c| std::cmp::Reverse(freq[c])); + let mut len_of = [0u8; 17]; + for (rank, &c) in order.iter().enumerate() { + len_of[c] = LENGTHS[rank]; + } + let mut counts = [0u8; 16]; + let mut vals: Vec = Vec::new(); + let mut code_of = [(0u32, 0u32); 17]; + let mut code = 0u32; + for len in 1..=16u8 { + for c in 0..17 { + if len_of[c] == len { + counts[len as usize - 1] += 1; + vals.push(c as u8); + code_of[c] = (code, u32::from(len)); + code += 1; + } + } + code <<= 1; + } + let mut out = vec![0xFF, 0xD8]; + // SOF3 + out.extend_from_slice(&[0xFF, 0xC3]); + out.extend_from_slice(&((8 + 3 * components) as u16).to_be_bytes()); + out.push(precision); + out.extend_from_slice(&(height as u16).to_be_bytes()); + out.extend_from_slice(&(width as u16).to_be_bytes()); + out.push(components as u8); + for c in 0..components { + out.extend_from_slice(&[c as u8 + 1, 0x11, 0]); + } + // DHT (one table, id 0) + out.extend_from_slice(&[0xFF, 0xC4]); + out.extend_from_slice(&((2 + 17 + vals.len()) as u16).to_be_bytes()); + out.push(0x00); + out.extend_from_slice(&counts); + out.extend_from_slice(&vals); + if let Some(r) = restart_rows { + out.extend_from_slice(&[0xFF, 0xDD, 0, 4]); + out.extend_from_slice(&((r * width) as u16).to_be_bytes()); + } + // SOS + out.extend_from_slice(&[0xFF, 0xDA]); + out.extend_from_slice(&((6 + 2 * components) as u16).to_be_bytes()); + out.push(components as u8); + for c in 0..components { + out.extend_from_slice(&[c as u8 + 1, 0x00]); + } + out.extend_from_slice(&[predictor, 0, 0]); + let mut w = BitWriter { out, acc: 0, n: 0 }; + let mut rst = 0u8; + for (i, &d) in diffs.iter().enumerate() { + if i > 0 && i % (interval.saturating_mul(row)) == 0 { + w.flush(); + w.out.extend_from_slice(&[0xFF, 0xD0 + rst]); + rst = (rst + 1) % 8; + } + let s = category(d); + let (c, l) = code_of[s as usize]; + w.put(c, l); + if (1..16).contains(&s) { + let extra = if d < 0 { (d - 1) as u32 & ((1 << s) - 1) } else { d as u32 }; + w.put(extra, s); + } + } + w.flush(); + let mut out = w.out; + out.extend_from_slice(&[0xFF, 0xD9]); + out +} + +// ---------------------------------------------------------------- TIFF writer + +/// A TIFF entry value. +#[derive(Debug, Clone)] +pub enum Val { + Byte(Vec), + Ascii(String), + Short(Vec), + Long(Vec), + Rational(Vec<(u32, u32)>), + SRational(Vec<(i32, i32)>), + Undefined(Vec), + /// Offsets of these blobs (LONG). + Blobs(Vec), + /// Offsets of these IFDs (LONG), for SubIFDs / EXIF. + Ifds(Vec), + /// An UNDEFINED value whose bytes are this IFD (a maker note). + IfdBytes(usize), +} + +/// A TIFF file under construction. +#[derive(Debug, Clone, Default)] +pub struct TiffBuilder { + pub big_endian: bool, + /// Entries per IFD (sorted by tag on write). + pub ifds: Vec>, + /// The main IFD chain, by index into `ifds`. + pub chain: Vec, + pub blobs: Vec>, + /// Write a CR2 header pointing at this IFD as the raw IFD. + pub cr2_raw_ifd: Option, +} + +impl TiffBuilder { + pub fn ifd(&mut self, entries: Vec<(u16, Val)>) -> usize { + self.ifds.push(entries); + self.ifds.len() - 1 + } + + pub fn blob(&mut self, data: Vec) -> usize { + self.blobs.push(data); + self.blobs.len() - 1 + } + + fn val_size(&self, v: &Val) -> usize { + match v { + Val::Byte(b) | Val::Undefined(b) => b.len(), + Val::Ascii(s) => s.len() + 1, + Val::Short(s) => 2 * s.len(), + Val::Long(l) => 4 * l.len(), + Val::Rational(r) => 8 * r.len(), + Val::SRational(r) => 8 * r.len(), + Val::Blobs(b) => 4 * b.len(), + Val::Ifds(i) => 4 * i.len(), + Val::IfdBytes(i) => self.ifd_size(*i), + } + } + + fn ifd_size(&self, i: usize) -> usize { + 2 + 12 * self.ifds[i].len() + 4 + } + + pub fn build(&self) -> Vec { + let header = if self.cr2_raw_ifd.is_some() { 16 } else { 8 }; + // Layout: each IFD, then its overflow values; then the blobs. + let mut pos = header; + let mut ifd_pos = vec![0usize; self.ifds.len()]; + let mut val_pos: Vec> = vec![Vec::new(); self.ifds.len()]; + // Maker-note IFDs are written inside their owner's value area. + let embedded: Vec = self.ifds.iter().flatten().filter_map(|(_, v)| if let Val::IfdBytes(i) = v { Some(*i) } else { None }).collect(); + for i in 0..self.ifds.len() { + if embedded.contains(&i) { + continue; + } + ifd_pos[i] = pos; + pos += self.ifd_size(i); + let mut vp = Vec::new(); + for (_, v) in self.sorted(i) { + let size = self.val_size(&v); + if size > 4 { + vp.push(pos); + if let Val::IfdBytes(e) = v { + ifd_pos[e] = pos; + // The embedded IFD's own values follow it. + let inner = pos + self.ifd_size(e); + let mut ip = Vec::new(); + let mut q = inner; + for (_, iv) in self.sorted(e) { + let s = self.val_size(&iv); + if s > 4 { + ip.push(q); + q += s + (s & 1); + } else { + ip.push(0); + } + } + val_pos[e] = ip; + pos = q; + } else { + pos += size + (size & 1); + } + } else { + vp.push(0); + } + } + val_pos[i] = vp; + } + let mut blob_pos = Vec::new(); + for b in &self.blobs { + blob_pos.push(pos); + pos += b.len() + (b.len() & 1); + } + let mut out = vec![0u8; pos]; + let be = self.big_endian; + let w16 = |o: &mut [u8], at: usize, v: u16| o[at..at + 2].copy_from_slice(&if be { v.to_be_bytes() } else { v.to_le_bytes() }); + let w32 = |o: &mut [u8], at: usize, v: u32| o[at..at + 4].copy_from_slice(&if be { v.to_be_bytes() } else { v.to_le_bytes() }); + out[0..4].copy_from_slice(if be { b"MM\0*" } else { b"II*\0" }); + w32(&mut out, 4, ifd_pos[self.chain[0]] as u32); + if let Some(r) = self.cr2_raw_ifd { + out[8..12].copy_from_slice(b"CR\x02\0"); + w32(&mut out, 12, ifd_pos[r] as u32); + } + for i in 0..self.ifds.len() { + let at = ifd_pos[i]; + let entries = self.sorted(i); + w16(&mut out, at, entries.len() as u16); + for (k, (tag, v)) in entries.iter().enumerate() { + let e = at + 2 + 12 * k; + let (typ, count, bytes) = self.encode(v, &blob_pos, &ifd_pos); + w16(&mut out, e, *tag); + w16(&mut out, e + 2, typ); + w32(&mut out, e + 4, count); + if bytes.len() <= 4 { + out[e + 8..e + 8 + bytes.len()].copy_from_slice(&bytes); + } else { + let vp = val_pos[i][k]; + w32(&mut out, e + 8, vp as u32); + if !matches!(v, Val::IfdBytes(_)) { + out[vp..vp + bytes.len()].copy_from_slice(&bytes); + } + } + } + let next = self.chain.iter().position(|&c| c == i).and_then(|p| self.chain.get(p + 1)).map(|&n| ifd_pos[n]).unwrap_or(0); + w32(&mut out, at + 2 + 12 * entries.len(), next as u32); + } + for (b, &p) in self.blobs.iter().zip(&blob_pos) { + out[p..p + b.len()].copy_from_slice(b); + } + out + } + + fn sorted(&self, i: usize) -> Vec<(u16, Val)> { + let mut e = self.ifds[i].clone(); + e.sort_by_key(|(t, _)| *t); + e + } + + fn encode(&self, v: &Val, blob_pos: &[usize], ifd_pos: &[usize]) -> (u16, u32, Vec) { + let be = self.big_endian; + let b16 = |v: u16| if be { v.to_be_bytes() } else { v.to_le_bytes() }; + let b32 = |v: u32| if be { v.to_be_bytes() } else { v.to_le_bytes() }; + match v { + Val::Byte(b) => (1, b.len() as u32, b.clone()), + Val::Undefined(b) => (7, b.len() as u32, b.clone()), + Val::Ascii(s) => { + let mut b = s.as_bytes().to_vec(); + b.push(0); + (2, b.len() as u32, b) + } + Val::Short(s) => (3, s.len() as u32, s.iter().flat_map(|&v| b16(v)).collect()), + Val::Long(l) => (4, l.len() as u32, l.iter().flat_map(|&v| b32(v)).collect()), + Val::Rational(r) => (5, r.len() as u32, r.iter().flat_map(|&(n, d)| [b32(n), b32(d)].concat()).collect()), + Val::SRational(r) => (10, r.len() as u32, r.iter().flat_map(|&(n, d)| [b32(n as u32), b32(d as u32)].concat()).collect()), + Val::Blobs(b) => (4, b.len() as u32, b.iter().flat_map(|&i| b32(blob_pos[i] as u32)).collect()), + Val::Ifds(i) => (4, i.len() as u32, i.iter().flat_map(|&k| b32(ifd_pos[k] as u32)).collect()), + Val::IfdBytes(i) => (7, self.ifd_size(*i) as u32, vec![0; 8]), + } + } +} + +fn srat(v: f64) -> (i32, i32) { + ((v * 10000.0).round() as i32, 10000) +} + +fn urat(v: f64) -> (u32, u32) { + ((v * 10000.0).round() as u32, 10000) +} + +// ---------------------------------------------------------------- DNG + +/// How the DNG stores its raw data. +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum DngStorage { + /// One strip per `rows` rows, uncompressed at `bits` per sample (8, 12, 16…). + Strips { rows: usize }, + /// Lossless-JPEG tiles of this size (2 components per JPEG row, as Adobe writes CFA tiles). + Lj92Tiles { width: usize, height: usize }, + /// Lossless-JPEG strips of this many rows, 1 component. + Lj92Strips { rows: usize }, +} + +/// A synthetic DNG. +#[derive(Debug, Clone)] +pub struct DngSpec { + pub width: usize, + pub height: usize, + /// 1 (CFA) or 3 (LinearRaw). + pub samples: usize, + /// Row-major samples (`width * height * samples`). + pub data: Vec, + pub bits: u16, + pub cfa: [u8; 4], + pub storage: DngStorage, + pub big_endian: bool, + pub black: Vec, + pub black_repeat: (u16, u16), + pub white: u32, + /// top, left, bottom, right + pub active_area: Option<[u32; 4]>, + /// origin (x, y), size (w, h), relative to the active area + pub default_crop: Option<([u32; 2], [u32; 2])>, + pub linearization: Option>, + /// (illuminant code, ColorMatrix row-major) + pub color_matrix1: Option<(u16, [f64; 9])>, + pub color_matrix2: Option<(u16, [f64; 9])>, + pub forward_matrix1: Option<[f64; 9]>, + pub forward_matrix2: Option<[f64; 9]>, + pub as_shot_neutral: Option<[f64; 3]>, + pub baseline_exposure: Option, + /// Raw OpcodeList2 bytes (see [`gain_map_opcode_list`]). + pub opcode_list2: Option>, + pub orientation: u16, + pub make: String, + pub model: String, +} + +impl DngSpec { + /// A CFA DNG with default metadata around `data` (RGGB, 16-bit, one strip). + pub fn cfa(width: usize, height: usize, data: Vec) -> Self { + DngSpec { + width, + height, + samples: 1, + data, + bits: 16, + cfa: [0, 1, 1, 2], + storage: DngStorage::Strips { rows: height }, + big_endian: false, + black: vec![0], + black_repeat: (1, 1), + white: 65535, + active_area: None, + default_crop: None, + linearization: None, + color_matrix1: None, + color_matrix2: None, + forward_matrix1: None, + forward_matrix2: None, + as_shot_neutral: None, + baseline_exposure: None, + orientation: 1, + opcode_list2: None, + make: "Photocraft".into(), + model: "Synthetic".into(), + } + } + + fn pack(&self, rows: std::ops::Range) -> Vec { + let n = self.width * self.samples; + let vals = &self.data[rows.start * n..rows.end * n]; + match self.bits { + 8 => vals.iter().map(|&v| v as u8).collect(), + 16 => vals.iter().flat_map(|&v| if self.big_endian { v.to_be_bytes() } else { v.to_le_bytes() }).collect(), + b => { + // MSB-first bit packing, each row padded to a byte. + let mut out = Vec::new(); + for row in vals.chunks(n) { + let mut acc: u64 = 0; + let mut have = 0u32; + for &v in row { + acc = (acc << b) | u64::from(v); + have += u32::from(b); + while have >= 8 { + have -= 8; + out.push((acc >> have) as u8); + } + } + if have > 0 { + out.push((acc << (8 - have)) as u8); + } + } + out + } + } + } + + pub fn build(&self) -> Vec { + let mut t = TiffBuilder { big_endian: self.big_endian, ..Default::default() }; + let mut raw: Vec<(u16, Val)> = vec![ + (254, Val::Long(vec![0])), + (256, Val::Long(vec![self.width as u32])), + (257, Val::Long(vec![self.height as u32])), + (258, Val::Short(vec![self.bits; self.samples])), + (262, Val::Short(vec![if self.samples == 1 { 32803 } else { 34892 }])), + (277, Val::Short(vec![self.samples as u16])), + (284, Val::Short(vec![1])), + (50717, Val::Long(vec![self.white; self.samples])), + (50713, Val::Short(vec![self.black_repeat.0, self.black_repeat.1])), + (50714, Val::Long(self.black.clone())), + ]; + if self.samples == 1 { + raw.push((33421, Val::Short(vec![2, 2]))); + raw.push((33422, Val::Byte(self.cfa.to_vec()))); + raw.push((50710, Val::Byte(vec![0, 1, 2]))); + raw.push((50711, Val::Short(vec![1]))); + } + match self.storage { + DngStorage::Strips { rows } => { + let mut offs = Vec::new(); + let mut lens = Vec::new(); + for y in (0..self.height).step_by(rows) { + let b = self.pack(y..(y + rows).min(self.height)); + lens.push(b.len() as u32); + offs.push(t.blob(b)); + } + raw.push((259, Val::Short(vec![1]))); + raw.push((278, Val::Long(vec![rows as u32]))); + raw.push((273, Val::Blobs(offs))); + raw.push((279, Val::Long(lens))); + } + DngStorage::Lj92Strips { rows } => { + let mut offs = Vec::new(); + let mut lens = Vec::new(); + for y in (0..self.height).step_by(rows) { + let r = (y + rows).min(self.height) - y; + let n = self.width * self.samples; + let j = lj92_encode(&self.data[y * n..(y + r) * n], self.width, r, self.samples, self.bits as u8, 1, None); + lens.push(j.len() as u32); + offs.push(t.blob(j)); + } + raw.push((259, Val::Short(vec![7]))); + raw.push((278, Val::Long(vec![rows as u32]))); + raw.push((273, Val::Blobs(offs))); + raw.push((279, Val::Long(lens))); + } + DngStorage::Lj92Tiles { width: tw, height: th } => { + let mut offs = Vec::new(); + let mut lens = Vec::new(); + let s = self.samples; + for ty in (0..self.height).step_by(th) { + for tx in (0..self.width).step_by(tw) { + // Edge tiles are padded by repeating the last row / column. + let mut tile = vec![0u16; tw * th * s]; + for y in 0..th { + for x in 0..tw { + let sy = (ty + y).min(self.height - 1); + let sx = (tx + x).min(self.width - 1); + for c in 0..s { + tile[(y * tw + x) * s + c] = self.data[(sy * self.width + sx) * s + c]; + } + } + } + // CFA tiles as Adobe writes them: half width, 2 components. + let j = if s == 1 && tw % 2 == 0 { + lj92_encode(&tile, tw / 2, th, 2, self.bits as u8, 1, None) + } else { + lj92_encode(&tile, tw, th, s, self.bits as u8, 1, None) + }; + lens.push(j.len() as u32); + offs.push(t.blob(j)); + } + } + raw.push((259, Val::Short(vec![7]))); + raw.push((322, Val::Long(vec![tw as u32]))); + raw.push((323, Val::Long(vec![th as u32]))); + raw.push((324, Val::Blobs(offs))); + raw.push((325, Val::Long(lens))); + } + } + if let Some(a) = self.active_area { + raw.push((50829, Val::Long(a.to_vec()))); + } + if let Some((o, s)) = self.default_crop { + raw.push((50719, Val::Long(o.to_vec()))); + raw.push((50720, Val::Long(s.to_vec()))); + } + if let Some(o) = &self.opcode_list2 { + raw.push((51009, Val::Undefined(o.clone()))); + } + if let Some(l) = &self.linearization { + raw.push((50712, Val::Short(l.clone()))); + } + let raw_ifd = t.ifd(raw); + // IFD0: a tiny RGB thumbnail plus the DNG metadata. + let thumb = t.blob(vec![128; 3 * 4]); + let mut ifd0: Vec<(u16, Val)> = vec![ + (254, Val::Long(vec![1])), + (256, Val::Long(vec![2])), + (257, Val::Long(vec![2])), + (258, Val::Short(vec![8, 8, 8])), + (259, Val::Short(vec![1])), + (262, Val::Short(vec![2])), + (271, Val::Ascii(self.make.clone())), + (272, Val::Ascii(self.model.clone())), + (273, Val::Blobs(vec![thumb])), + (274, Val::Short(vec![self.orientation])), + (277, Val::Short(vec![3])), + (278, Val::Long(vec![2])), + (279, Val::Long(vec![12])), + (330, Val::Ifds(vec![raw_ifd])), + (50706, Val::Byte(vec![1, 4, 0, 0])), + (50708, Val::Ascii(format!("{} {}", self.make, self.model))), + ]; + let m = |v: &[f64; 9]| Val::SRational(v.iter().map(|&x| srat(x)).collect()); + if let Some((ill, cm)) = &self.color_matrix1 { + ifd0.push((50721, m(cm))); + ifd0.push((50778, Val::Short(vec![*ill]))); + } + if let Some((ill, cm)) = &self.color_matrix2 { + ifd0.push((50722, m(cm))); + ifd0.push((50779, Val::Short(vec![*ill]))); + } + if let Some(f) = &self.forward_matrix1 { + ifd0.push((50964, m(f))); + } + if let Some(f) = &self.forward_matrix2 { + ifd0.push((50965, m(f))); + } + if let Some(n) = self.as_shot_neutral { + ifd0.push((50728, Val::Rational(n.iter().map(|&v| urat(v)).collect()))); + } + if let Some(b) = self.baseline_exposure { + ifd0.push((50730, Val::SRational(vec![srat(b)]))); + } + let i0 = t.ifd(ifd0); + t.chain = vec![i0]; + t.build() + } +} + +/// An OpcodeList2 with one GainMap over active-area rectangle `area` (top, +/// left, bottom, right), every `pitch` rows and columns, from a 2×2 grid of +/// gains spanning the image (relative coordinates 0..1). +pub fn gain_map_opcode_list(area: [u32; 4], pitch: u32, gains: [[f32; 2]; 2]) -> Vec { + let mut p = Vec::new(); + for v in [area[0], area[1], area[2], area[3], 0, 1, pitch, pitch, 2, 2] { + p.extend_from_slice(&v.to_be_bytes()); + } + for v in [1.0f64, 1.0, 0.0, 0.0] { + p.extend_from_slice(&v.to_be_bytes()); + } + p.extend_from_slice(&1u32.to_be_bytes()); + for g in gains.iter().flatten() { + p.extend_from_slice(&g.to_bits().to_be_bytes()); + } + let mut b = 1u32.to_be_bytes().to_vec(); + for v in [9u32, 0x0103_0000, 1, p.len() as u32] { + b.extend_from_slice(&v.to_be_bytes()); + } + b.extend_from_slice(&p); + b +} + +// ---------------------------------------------------------------- CR2 + +/// A synthetic CR2: one lossless JPEG with `components` components, cut into +/// vertical slices, plus a Canon maker note with SensorInfo and ColorData. +#[derive(Debug, Clone)] +pub struct Cr2Spec { + pub width: usize, + pub height: usize, + /// Sensor samples, row-major, RGGB from (`left`, `top`). + pub data: Vec, + pub precision: u8, + pub components: usize, + /// Slice widths (all but the last equal); empty = no slice tag. + pub slices: Vec, + /// Image borders: left, top, right, bottom (inclusive), as in SensorInfo. + pub borders: Option<[u16; 4]>, + /// As-shot RGGB levels written to ColorData at word offset 0x3F. + pub wb_rggb: Option<[u16; 4]>, + pub orientation: u16, +} + +impl Cr2Spec { + pub fn build(&self) -> Vec { + let (w, h) = (self.width, self.height); + // Slice the sensor data into the JPEG sample order. + let widths: Vec = if self.slices.is_empty() { vec![w] } else { self.slices.clone() }; + assert_eq!(widths.iter().sum::(), w); + let mut seq = Vec::with_capacity(w * h); + let mut x0 = 0; + for &sw in &widths { + for y in 0..h { + seq.extend_from_slice(&self.data[y * w + x0..y * w + x0 + sw]); + } + x0 += sw; + } + let c = self.components; + assert_eq!(w % c, 0); + let jpeg = lj92_encode(&seq, w / c, h, c, self.precision, 1, None); + let mut t = TiffBuilder::default(); + let jb = t.blob(jpeg.clone()); + let mut raw = vec![(259, Val::Short(vec![6])), (273, Val::Blobs(vec![jb])), (279, Val::Long(vec![jpeg.len() as u32]))]; + if !self.slices.is_empty() { + let n = widths.len() - 1; + raw.push((50752, Val::Short(vec![n as u16, widths[0] as u16, widths[n] as u16]))); + } + let raw_ifd = t.ifd(raw); + let mut mn: Vec<(u16, Val)> = Vec::new(); + if let Some([l, tp, r, b]) = self.borders { + mn.push((0x00E0, Val::Short(vec![34, w as u16, h as u16, 0, 0, l, tp, r, b, 0, 0, 0, 0, 0, 0, 0, 0]))); + } + if let Some(wb) = self.wb_rggb { + let mut cd = vec![0u16; 1273]; + cd[0x3F..0x43].copy_from_slice(&wb); + mn.push((0x4001, Val::Short(cd))); + } + let mut exif: Vec<(u16, Val)> = vec![(33434, Val::Rational(vec![(1, 100)]))]; + if !mn.is_empty() { + let mi = t.ifd(mn); + exif.push((37500, Val::IfdBytes(mi))); + } + let exif_ifd = t.ifd(exif); + let thumb = t.blob(vec![128; 12]); + let ifd0 = t.ifd(vec![ + (256, Val::Long(vec![2])), + (257, Val::Long(vec![2])), + (259, Val::Short(vec![1])), + (262, Val::Short(vec![2])), + (271, Val::Ascii("Canon".into())), + (272, Val::Ascii("Canon EOS Synthetic".into())), + (273, Val::Blobs(vec![thumb])), + (274, Val::Short(vec![self.orientation])), + (277, Val::Short(vec![3])), + (279, Val::Long(vec![12])), + (34665, Val::Ifds(vec![exif_ifd])), + ]); + let ifd1 = t.ifd(vec![(256, Val::Long(vec![1]))]); + let ifd2 = t.ifd(vec![(256, Val::Long(vec![1]))]); + t.chain = vec![ifd0, ifd1, ifd2, raw_ifd]; + t.cr2_raw_ifd = Some(raw_ifd); + t.build() + } +} + +// ---------------------------------------------------------------- TIFF/EP + +/// A synthetic uncompressed TIFF/EP raw (NEF / ARW-like): IFD0 with Make, the +/// CFA image in a SubIFD with 16-bit samples. +pub fn tiff_ep(make: &str, width: usize, height: usize, data: &[u16], cfa: [u8; 4], bits: u16, extra_raw_tags: Vec<(u16, Val)>) -> Vec { + let mut t = TiffBuilder::default(); + let strip = t.blob(data.iter().flat_map(|v| v.to_le_bytes()).collect()); + let mut raw = vec![ + (254, Val::Long(vec![0])), + (256, Val::Long(vec![width as u32])), + (257, Val::Long(vec![height as u32])), + (258, Val::Short(vec![bits])), + (259, Val::Short(vec![1])), + (262, Val::Short(vec![32803])), + (273, Val::Blobs(vec![strip])), + (277, Val::Short(vec![1])), + (278, Val::Long(vec![height as u32])), + (279, Val::Long(vec![(data.len() * 2) as u32])), + (33421, Val::Short(vec![2, 2])), + (33422, Val::Byte(cfa.to_vec())), + ]; + raw.extend(extra_raw_tags); + let raw_ifd = t.ifd(raw); + let thumb = t.blob(vec![128; 12]); + let ifd0 = t.ifd(vec![ + (254, Val::Long(vec![1])), + (256, Val::Long(vec![2])), + (257, Val::Long(vec![2])), + (259, Val::Short(vec![1])), + (262, Val::Short(vec![2])), + (271, Val::Ascii(make.into())), + (273, Val::Blobs(vec![thumb])), + (277, Val::Short(vec![3])), + (279, Val::Long(vec![12])), + (330, Val::Ifds(vec![raw_ifd])), + ]); + t.chain = vec![ifd0]; + t.build() +} + +// ---------------------------------------------------------------- scenes + +/// A smooth, colourful synthetic scene as linear RGB in 0..1. +pub fn scene(width: usize, height: usize) -> Vec<[f32; 3]> { + (0..width * height) + .map(|i| { + let (x, y) = ((i % width) as f32 / width.max(1) as f32, (i / width) as f32 / height.max(1) as f32); + [0.1 + 0.7 * x, 0.15 + 0.5 * (1.0 - y) * x + 0.2 * y, 0.6 - 0.4 * y + 0.2 * x] + }) + .collect() +} + +/// Mosaics linear RGB through a 2×2 CFA (`cfa` row-major colours) into +/// sensor values `black + v * (white - black)`. +pub fn mosaic(rgb: &[[f32; 3]], width: usize, cfa: [u8; 4], black: u16, white: u16) -> Vec { + rgb.iter() + .enumerate() + .map(|(i, p)| { + let (x, y) = (i % width, i / width); + let c = cfa[(y & 1) * 2 + (x & 1)] as usize; + let v = p[c].clamp(0.0, 1.0); + (f32::from(black) + v * f32::from(white - black)).round() as u16 + }) + .collect() +} diff --git a/crates/raw/src/tiff.rs b/crates/raw/src/tiff.rs new file mode 100644 index 00000000..764cc870 --- /dev/null +++ b/crates/raw/src/tiff.rs @@ -0,0 +1,381 @@ +//! A small, bounds-checked TIFF 6.0 / TIFF-EP structure reader. +//! +//! Only what raw decoding needs: IFD entries with typed value access, the IFD +//! chain, SubIFDs (tag 330), the EXIF IFD (34665) and IFDs embedded at other +//! offsets (maker notes). Every offset is checked against the buffer; reads +//! that fall outside return `None`. The walk is bounded in IFD count and +//! entry count so hostile files cannot loop or blow up memory. + +/// Tag numbers (TIFF 6.0, TIFF/EP, EXIF and DNG 1.7). +pub(crate) mod tag { + pub const NEW_SUBFILE_TYPE: u16 = 254; + pub const IMAGE_WIDTH: u16 = 256; + pub const IMAGE_LENGTH: u16 = 257; + pub const BITS_PER_SAMPLE: u16 = 258; + pub const COMPRESSION: u16 = 259; + pub const PHOTOMETRIC: u16 = 262; + pub const MAKE: u16 = 271; + pub const MODEL: u16 = 272; + pub const STRIP_OFFSETS: u16 = 273; + pub const ORIENTATION: u16 = 274; + pub const SAMPLES_PER_PIXEL: u16 = 277; + pub const ROWS_PER_STRIP: u16 = 278; + pub const STRIP_BYTE_COUNTS: u16 = 279; + pub const PLANAR_CONFIGURATION: u16 = 284; + pub const TILE_WIDTH: u16 = 322; + pub const TILE_LENGTH: u16 = 323; + pub const TILE_OFFSETS: u16 = 324; + pub const TILE_BYTE_COUNTS: u16 = 325; + pub const SUB_IFDS: u16 = 330; + pub const SAMPLE_FORMAT: u16 = 339; + pub const JPEG_INTERCHANGE_FORMAT: u16 = 513; + pub const JPEG_INTERCHANGE_FORMAT_LENGTH: u16 = 514; + pub const CFA_REPEAT_PATTERN_DIM: u16 = 33421; + pub const CFA_PATTERN: u16 = 33422; + pub const EXIF_IFD: u16 = 34665; + pub const EXIF_CFA_PATTERN: u16 = 41730; + pub const MAKER_NOTE: u16 = 37500; + pub const DNG_VERSION: u16 = 50706; + pub const UNIQUE_CAMERA_MODEL: u16 = 50708; + pub const CFA_PLANE_COLOR: u16 = 50710; + pub const CFA_LAYOUT: u16 = 50711; + pub const LINEARIZATION_TABLE: u16 = 50712; + pub const BLACK_LEVEL_REPEAT_DIM: u16 = 50713; + pub const BLACK_LEVEL: u16 = 50714; + pub const BLACK_LEVEL_DELTA_H: u16 = 50715; + pub const BLACK_LEVEL_DELTA_V: u16 = 50716; + pub const WHITE_LEVEL: u16 = 50717; + pub const DEFAULT_SCALE: u16 = 50718; + pub const DEFAULT_CROP_ORIGIN: u16 = 50719; + pub const DEFAULT_CROP_SIZE: u16 = 50720; + pub const COLOR_MATRIX_1: u16 = 50721; + pub const COLOR_MATRIX_2: u16 = 50722; + pub const CAMERA_CALIBRATION_1: u16 = 50723; + pub const CAMERA_CALIBRATION_2: u16 = 50724; + pub const ANALOG_BALANCE: u16 = 50727; + pub const AS_SHOT_NEUTRAL: u16 = 50728; + pub const AS_SHOT_WHITE_XY: u16 = 50729; + pub const BASELINE_EXPOSURE: u16 = 50730; + pub const CALIBRATION_ILLUMINANT_1: u16 = 50778; + pub const CALIBRATION_ILLUMINANT_2: u16 = 50779; + pub const ACTIVE_AREA: u16 = 50829; + pub const FORWARD_MATRIX_1: u16 = 50964; + pub const FORWARD_MATRIX_2: u16 = 50965; + pub const OPCODE_LIST_1: u16 = 51008; + pub const OPCODE_LIST_2: u16 = 51009; + pub const OPCODE_LIST_3: u16 = 51022; + pub const BASELINE_EXPOSURE_OFFSET: u16 = 51109; + // CR2 private tag in the raw IFD: slice count, slice width, last slice width. + pub const CR2_SLICE: u16 = 50752; +} + +/// Most IFDs a walk will visit. +const MAX_IFDS: usize = 64; +/// Most entries read from one IFD. +const MAX_ENTRIES: usize = 4096; +/// Most values decoded from one entry by the list accessors. +pub(crate) const MAX_VALUES: usize = 1 << 20; + +/// One directory entry. +#[derive(Debug, Clone, Copy)] +pub(crate) struct Entry { + pub tag: u16, + pub typ: u16, + pub count: u32, + /// Absolute offset of the value bytes in the buffer (may be out of range; + /// accessors check). + pub at: usize, +} + +impl Entry { + /// Size in bytes of one value of this entry's type (0 = unknown type). + fn unit(&self) -> usize { + type_size(self.typ) + } +} + +/// One image file directory. +#[derive(Debug, Clone, Default)] +pub(crate) struct Ifd { + pub offset: usize, + pub entries: Vec, +} + +impl Ifd { + pub fn get(&self, tag: u16) -> Option<&Entry> { + self.entries.iter().find(|e| e.tag == tag) + } + pub fn has(&self, tag: u16) -> bool { + self.get(tag).is_some() + } +} + +pub(crate) fn type_size(typ: u16) -> usize { + match typ { + 1 | 2 | 6 | 7 => 1, + 3 | 8 => 2, + 4 | 9 | 11 | 13 => 4, + 5 | 10 | 12 => 8, + _ => 0, + } +} + +/// A TIFF-structured buffer. +#[derive(Debug, Clone, Copy)] +pub(crate) struct Tiff<'a> { + pub data: &'a [u8], + pub le: bool, + /// Offset of IFD0 from the header. + pub first_ifd: usize, +} + +impl<'a> Tiff<'a> { + /// Parses the 8-byte header. Accepts the TIFF magics and the TIFF-like + /// variants used by camera raws (ORF `IIRO`/`IIRS`/`MMOR`, RW2 `IIU\0`). + pub fn new(data: &'a [u8]) -> Option { + let head = data.get(0..4)?; + let le = match head { + b"II*\0" | b"IIRO" | b"IIRS" | b"IIU\0" => true, + b"MM\0*" | b"MMOR" => false, + _ => return None, + }; + let mut t = Tiff { data, le, first_ifd: 0 }; + t.first_ifd = t.u32_at(4)? as usize; + Some(t) + } + + pub fn bytes(&self, at: usize, len: usize) -> Option<&'a [u8]> { + self.data.get(at..at.checked_add(len)?) + } + + pub fn u16_at(&self, at: usize) -> Option { + let s: [u8; 2] = self.bytes(at, 2)?.try_into().ok()?; + Some(if self.le { u16::from_le_bytes(s) } else { u16::from_be_bytes(s) }) + } + + pub fn u32_at(&self, at: usize) -> Option { + let s: [u8; 4] = self.bytes(at, 4)?.try_into().ok()?; + Some(if self.le { u32::from_le_bytes(s) } else { u32::from_be_bytes(s) }) + } + + fn u64_at(&self, at: usize) -> Option { + let s: [u8; 8] = self.bytes(at, 8)?.try_into().ok()?; + Some(if self.le { u64::from_le_bytes(s) } else { u64::from_be_bytes(s) }) + } + + /// Reads the IFD at `offset`. Value offsets are relative to `base` + /// (0 for ordinary IFDs; maker notes sometimes use their own base). + pub fn ifd_at(&self, offset: usize, base: usize) -> Option { + let n = usize::from(self.u16_at(offset)?); + if n == 0 || n > MAX_ENTRIES { + return None; + } + let mut entries = Vec::with_capacity(n); + for i in 0..n { + let e = offset.checked_add(2 + 12 * i)?; + let (Some(tag), Some(typ), Some(count)) = (self.u16_at(e), self.u16_at(e.saturating_add(2)), self.u32_at(e.saturating_add(4))) else { + break; + }; + let size = type_size(typ).saturating_mul(count as usize); + let at = if size <= 4 { e.saturating_add(8) } else { (self.u32_at(e.saturating_add(8))? as usize).saturating_add(base) }; + entries.push(Entry { tag, typ, count, at }); + } + Some(Ifd { offset, entries }) + } + + /// Offset of the IFD that follows the one at `ifd` (0 = none). + pub fn next_ifd(&self, ifd: &Ifd) -> usize { + let n = ifd.entries.len(); + self.u32_at(ifd.offset.saturating_add(2 + 12 * n)).unwrap_or(0) as usize + } + + /// Every IFD reachable from the header: the main chain, SubIFDs and the + /// EXIF IFD, depth-first, deduplicated and bounded. + pub fn all_ifds(&self) -> Vec { + let mut out: Vec = Vec::new(); + let mut seen: Vec = Vec::new(); + let mut stack = vec![self.first_ifd]; + while let Some(off) = stack.pop() { + if off == 0 || seen.contains(&off) || seen.len() >= MAX_IFDS { + continue; + } + seen.push(off); + let Some(ifd) = self.ifd_at(off, 0) else { continue }; + let next = self.next_ifd(&ifd); + if next != 0 { + stack.push(next); + } + for t in [tag::SUB_IFDS, tag::EXIF_IFD] { + if let Some(e) = ifd.get(t) { + // Reverse so the first SubIFD is visited first. + let mut v = self.uints(e); + v.truncate(16); + stack.extend(v.into_iter().rev().map(|o| o as usize)); + } + } + out.push(ifd); + } + out + } + + /// Value `i` of an entry as an unsigned integer (BYTE/SHORT/LONG/IFD and + /// their signed variants, negative values are rejected). + pub fn uint(&self, e: &Entry, i: usize) -> Option { + if i >= e.count as usize { + return None; + } + let at = e.at.checked_add(i.checked_mul(e.unit())?)?; + match e.typ { + 1 | 7 => self.bytes(at, 1).map(|b| u32::from(b[0])), + 6 => self.bytes(at, 1).and_then(|b| u32::try_from(b[0] as i8).ok()), + 3 => self.u16_at(at).map(u32::from), + 8 => self.u16_at(at).and_then(|v| u32::try_from(v as i16).ok()), + 4 | 13 => self.u32_at(at), + 9 => self.u32_at(at).and_then(|v| u32::try_from(v as i32).ok()), + _ => None, + } + } + + /// Value `i` as a float (any numeric type, rationals included). + pub fn float(&self, e: &Entry, i: usize) -> Option { + if i >= e.count as usize { + return None; + } + let at = e.at.checked_add(i.checked_mul(e.unit())?)?; + let v = match e.typ { + 1 | 7 => f64::from(self.bytes(at, 1)?[0]), + 6 => f64::from(self.bytes(at, 1)?[0] as i8), + 3 => f64::from(self.u16_at(at)?), + 8 => f64::from(self.u16_at(at)? as i16), + 4 | 13 => f64::from(self.u32_at(at)?), + 9 => f64::from(self.u32_at(at)? as i32), + 5 => { + let (n, d) = (self.u32_at(at)?, self.u32_at(at.saturating_add(4))?); + if d == 0 { + return None; + } + f64::from(n) / f64::from(d) + } + 10 => { + let (n, d) = (self.u32_at(at)? as i32, self.u32_at(at.saturating_add(4))? as i32); + if d == 0 { + return None; + } + f64::from(n) / f64::from(d) + } + 11 => f64::from(f32::from_bits(self.u32_at(at)?)), + 12 => f64::from_bits(self.u64_at(at)?), + _ => return None, + }; + v.is_finite().then_some(v) + } + + /// All values as unsigned integers (empty if any is unreadable). + pub fn uints(&self, e: &Entry) -> Vec { + let n = (e.count as usize).min(MAX_VALUES); + (0..n).map(|i| self.uint(e, i)).collect::>>().unwrap_or_default() + } + + /// All values as floats (empty if any is unreadable). + pub fn floats(&self, e: &Entry) -> Vec { + let n = (e.count as usize).min(MAX_VALUES); + (0..n).map(|i| self.float(e, i)).collect::>>().unwrap_or_default() + } + + /// The raw value bytes of an entry. + pub fn raw(&self, e: &Entry) -> Option<&'a [u8]> { + let len = e.unit().checked_mul(e.count as usize)?; + self.bytes(e.at, len) + } + + /// An ASCII value, trimmed of NULs and whitespace. + pub fn ascii(&self, e: &Entry) -> Option { + let b = self.raw(e)?; + let end = b.iter().position(|&c| c == 0).unwrap_or(b.len()); + let s = String::from_utf8_lossy(&b[..end]).trim().to_string(); + (!s.is_empty()).then_some(s) + } + + pub fn tag_uint(&self, ifd: &Ifd, t: u16) -> Option { + ifd.get(t).and_then(|e| self.uint(e, 0)) + } + + pub fn tag_uints(&self, ifd: &Ifd, t: u16) -> Vec { + ifd.get(t).map(|e| self.uints(e)).unwrap_or_default() + } + + pub fn tag_floats(&self, ifd: &Ifd, t: u16) -> Vec { + ifd.get(t).map(|e| self.floats(e)).unwrap_or_default() + } + + pub fn tag_ascii(&self, ifd: &Ifd, t: u16) -> Option { + ifd.get(t).and_then(|e| self.ascii(e)) + } +} + +/// A human-readable listing of every IFD and entry (debugging aid). +pub(crate) fn dump(data: &[u8]) -> String { + use std::fmt::Write; + let Some(t) = Tiff::new(data) else { return "not a TIFF-structured file".into() }; + let mut s = String::new(); + for ifd in t.all_ifds() { + let _ = writeln!(s, "IFD @{}", ifd.offset); + for e in &ifd.entries { + let preview: String = if e.typ == 2 { + t.ascii(e).unwrap_or_default() + } else { + let n = (e.count as usize).min(12); + let v: Vec = (0..n).map(|i| t.float(e, i).map(|f| format!("{f}")).unwrap_or_else(|| "?".into())).collect(); + format!("{}{}", v.join(" "), if e.count as usize > n { " …" } else { "" }) + }; + let _ = writeln!(s, " {:5} (0x{:04X}) type {:2} count {:7} @{:9}: {}", e.tag, e.tag, e.typ, e.count, e.at, preview); + } + } + s +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn header_variants() { + assert!(Tiff::new(b"II*\0\x08\0\0\0").is_some()); + assert!(Tiff::new(b"MM\0*\0\0\0\x08").is_some()); + assert!(Tiff::new(b"IIRO\x08\0\0\0").is_some()); + assert!(Tiff::new(b"IIU\0\x08\0\0\0").is_some()); + assert!(Tiff::new(b"II*\0").is_none()); + assert!(Tiff::new(b"PNG\0\0\0\0\0").is_none()); + } + + #[test] + fn self_referencing_ifds_terminate() { + // IFD0 at 8 with one entry whose SubIFD points at itself, and a next + // pointer that also points at itself. + let mut b = b"II*\0\x08\0\0\0".to_vec(); + b.extend_from_slice(&1u16.to_le_bytes()); + b.extend_from_slice(&330u16.to_le_bytes()); + b.extend_from_slice(&4u16.to_le_bytes()); + b.extend_from_slice(&1u32.to_le_bytes()); + b.extend_from_slice(&8u32.to_le_bytes()); + b.extend_from_slice(&8u32.to_le_bytes()); + let t = Tiff::new(&b).unwrap(); + assert_eq!(t.all_ifds().len(), 1); + } + + #[test] + fn rational_with_zero_denominator_is_none() { + let mut b = b"II*\0\x08\0\0\0".to_vec(); + b.extend_from_slice(&1u16.to_le_bytes()); + b.extend_from_slice(&50730u16.to_le_bytes()); + b.extend_from_slice(&10u16.to_le_bytes()); + b.extend_from_slice(&1u32.to_le_bytes()); + b.extend_from_slice(&26u32.to_le_bytes()); + b.extend_from_slice(&0u32.to_le_bytes()); + b.extend_from_slice(&5i32.to_le_bytes()); + b.extend_from_slice(&0i32.to_le_bytes()); + let t = Tiff::new(&b).unwrap(); + let ifd = t.ifd_at(8, 0).unwrap(); + assert_eq!(t.float(ifd.get(50730).unwrap(), 0), None); + } +} diff --git a/crates/raw/src/tiffep.rs b/crates/raw/src/tiffep.rs new file mode 100644 index 00000000..d7f7be88 --- /dev/null +++ b/crates/raw/src/tiffep.rs @@ -0,0 +1,177 @@ +//! TIFF/EP-structured raws (NEF, ARW, PEF and others) whose sensor data is +//! stored uncompressed or as lossless JPEG in a standard CFA IFD (TIFF/EP, +//! ISO 12234-2: PhotometricInterpretation 32803, CFARepeatPatternDim and +//! CFAPattern). Vendor-specific compressions are reported as unsupported. +//! +//! Black and white levels use DNG-style tags when the camera writes them, +//! then the vendor's publicly documented tags (Nikon maker note BlackLevel, +//! Sony raw-IFD BlackLevel), else 0 with a warning. The white level falls +//! back to the clipping point found in the data. The as-shot white balance +//! comes from the documented vendor tags (Nikon WB_RBLevels, Sony +//! WB_RGGBLevels); otherwise it is estimated. + +use crate::cr2::clip_level; +use crate::error::{RawError, Result}; +use crate::sensor::{BlackLevels, Cfa, JpegLayout, Rect, Sensor, read_plane}; +use crate::tiff::{Ifd, Tiff, tag}; +use crate::{Limits, RawFormat}; + +const PHOTOMETRIC_CFA: u32 = 32803; +const NIKON_WB_RB_LEVELS: u16 = 0x000C; +const NIKON_BLACK_LEVEL: u16 = 0x003D; +const SONY_RAW_FILE_TYPE: u16 = 0x7000; +const SONY_BLACK_LEVEL: u16 = 0x7310; +const SONY_WB_RGGB_LEVELS: u16 = 0x7313; + +fn raw_ifd(t: &Tiff, ifds: &[Ifd]) -> Option { + ifds.iter() + .filter(|i| t.tag_uint(i, tag::PHOTOMETRIC) == Some(PHOTOMETRIC_CFA)) + .max_by_key(|i| u64::from(t.tag_uint(i, tag::IMAGE_WIDTH).unwrap_or(0)) * u64::from(t.tag_uint(i, tag::IMAGE_LENGTH).unwrap_or(0))) + .cloned() +} + +/// `true` for Sony ARW files whose raw IFD is not CFA data (the reduced-size +/// lossless "M" / "S" YCbCr variants): recognised, but not decoded. +pub(crate) fn is_sony_non_cfa(t: &Tiff, ifds: &[Ifd]) -> bool { + ifds.iter().any(|i| i.has(SONY_RAW_FILE_TYPE)) && ifds.iter().any(|i| t.tag_ascii(i, tag::MAKE).is_some_and(|m| m.to_ascii_uppercase().starts_with("SONY"))) +} + +/// CFA from the TIFF/EP tags of the raw IFD, or the EXIF CFAPattern. +fn cfa(t: &Tiff, raw: &Ifd, ifds: &[Ifd]) -> Option { + let dim = t.tag_uints(raw, tag::CFA_REPEAT_PATTERN_DIM); + let pat = t.tag_uints(raw, tag::CFA_PATTERN); + let (rows, cols, colors) = match (dim.as_slice(), pat) { + ([r, c], p) if p.len() == (*r as usize) * (*c as usize) && !p.is_empty() => (*r as usize, *c as usize, p), + _ => { + // EXIF CFAPattern: two 16-bit counts (columns, rows) then the colours. + let e = ifds.iter().find_map(|i| i.get(tag::EXIF_CFA_PATTERN).copied())?; + let b = t.raw(&e)?; + let rd = |o: usize| b.get(o..o + 2).map(|s| if t.le { u16::from_le_bytes([s[0], s[1]]) } else { u16::from_be_bytes([s[0], s[1]]) }); + let (c, r) = (usize::from(rd(0)?), usize::from(rd(2)?)); + let p = b.get(4..4 + r.checked_mul(c)?)?; + (r, c, p.iter().map(|&v| u32::from(v)).collect()) + } + }; + if rows == 0 || cols == 0 || rows > 16 || cols > 16 { + return None; + } + let colors = colors.into_iter().map(|c| u8::try_from(c).ok().filter(|c| *c <= 2)).collect::>>()?; + Some(Cfa { width: cols, height: rows, colors, origin_x: 0, origin_y: 0 }) +} + +/// Nikon maker note: "Nikon\0", version, then an embedded TIFF header whose +/// offsets are relative to itself. +fn nikon_maker_note<'a>(t: &Tiff<'a>, ifds: &[Ifd]) -> Option<(Tiff<'a>, Ifd)> { + let e = ifds.iter().find_map(|i| i.get(tag::MAKER_NOTE).copied())?; + let b = t.raw(&e)?; + if !b.starts_with(b"Nikon\0") { + return None; + } + let base = e.at.checked_add(10)?; + let inner = Tiff::new(t.data.get(base..)?)?; + let ifd = inner.ifd_at(inner.first_ifd, 0)?; + Some((inner, ifd)) +} + +/// Four levels listed in R, G, G, B order, placed at their 2×2 CFA positions +/// (data parity). +fn rggb_by_position(cfa: &Cfa, v: &[f64]) -> Option<[f32; 4]> { + let [r, g1, g2, b] = <[f64; 4]>::try_from(v).ok()?; + if [r, g1, g2, b].iter().any(|x| !x.is_finite() || *x < 0.0) { + return None; + } + let mut greens = [g1, g2].into_iter(); + let mut out = [0.0f32; 4]; + for (k, slot) in out.iter_mut().enumerate() { + *slot = match cfa.color(k & 1, k >> 1) { + 0 => r, + 2 => b, + _ => greens.next().unwrap_or(g1), + } as f32; + } + Some(out) +} + +pub(crate) fn decode(t: &Tiff, format: RawFormat, limits: &Limits) -> Result { + let ifds = t.all_ifds(); + let Some(raw) = raw_ifd(t, &ifds) else { + if is_sony_non_cfa(t, &ifds) { + return Err(RawError::unsupported("Sony reduced-size (M / S) lossless ARW")); + } + return Err(RawError::unsupported(format!("{}: no CFA image found", format.name()))); + }; + let plane = read_plane(t, &raw, limits, JpegLayout::Quads)?; + if plane.samples != 1 { + return Err(RawError::unsupported(format!("CFA data with {} samples per pixel", plane.samples))); + } + let cfa = cfa(t, &raw, &ifds).ok_or_else(|| RawError::unsupported(format!("{}: CFA pattern not recorded", format.name())))?; + if !cfa.is_bayer() { + return Err(RawError::unsupported(format!("{}x{} non-Bayer CFA", cfa.height, cfa.width))); + } + let mut warnings = Vec::new(); + let full = Rect::new(0, 0, plane.width, plane.height); + let make = ifds.iter().find_map(|i| t.tag_ascii(i, tag::MAKE)); + let model = ifds.iter().find_map(|i| t.tag_ascii(i, tag::MODEL)); + + let nikon = if format == RawFormat::Nef { nikon_maker_note(t, &ifds) } else { None }; + let bl = t.tag_floats(&raw, tag::BLACK_LEVEL); + let vendor_black = match &nikon { + Some((n, m)) => n.tag_floats(m, NIKON_BLACK_LEVEL), + None => t.tag_floats(&raw, SONY_BLACK_LEVEL), + }; + let black = if let Some(b) = bl.first() { + BlackLevels::uniform(*b as f32) + } else if let Some(v) = rggb_by_position(&cfa, &vendor_black) { + BlackLevels { rows: 2, cols: 2, values: v.to_vec(), delta_h: Vec::new(), delta_v: Vec::new() } + } else { + warnings.push(format!("{}: black level not recorded in a documented tag; assumed 0", format.name())); + BlackLevels::uniform(0.0) + }; + let wl = t.tag_floats(&raw, tag::WHITE_LEVEL); + let white = match wl.first() { + Some(w) if *w > 0.0 => *w as f32, + _ => clip_level(&plane.data, plane.bits), + }; + let camera_wb = match &nikon { + Some((n, m)) => match n.tag_floats(m, NIKON_WB_RB_LEVELS).as_slice() { + [r, b, ..] => Some([*r, 1.0, *b]), + _ => None, + }, + None => match t.tag_floats(&raw, SONY_WB_RGGB_LEVELS).as_slice() { + [r, g1, g2, b] if *g1 + *g2 > 0.0 => Some([2.0 * r / (g1 + g2), 1.0, 2.0 * b / (g1 + g2)]), + _ => None, + }, + } + .filter(|m| m.iter().all(|v| (0.25..8.0).contains(v))); + // DNG-style DefaultCropOrigin / DefaultCropSize (some vendors write them in the raw IFD). + let mut crop = full; + if let ([ox, oy], [cw, ch]) = (t.tag_floats(&raw, tag::DEFAULT_CROP_ORIGIN).as_slice(), t.tag_floats(&raw, tag::DEFAULT_CROP_SIZE).as_slice()) { + let f = |v: f64| if v.is_finite() && v >= 0.0 { v.round() as usize } else { 0 }; + let r = Rect::new(f(*ox), f(*oy), f(*cw), f(*ch)).intersect(&full); + if !r.is_empty() { + crop = r; + } + } + let ifd0 = t.ifd_at(t.first_ifd, 0); + Ok(Sensor { + format, + make, + model, + width: plane.width, + height: plane.height, + samples: 1, + data: plane.data, + cfa: Some(cfa), + linearization: None, + black, + white: [white; 3], + active: full, + crop, + color: Default::default(), + camera_wb, + orientation: ifd0.and_then(|i| t.tag_uint(&i, tag::ORIENTATION)).map(|o| o as u16).filter(|o| (1..=8).contains(o)).unwrap_or(1), + baseline_exposure: 0.0, + gain_maps: Vec::new(), + warnings, + }) +} diff --git a/crates/raw/tests/dng.rs b/crates/raw/tests/dng.rs new file mode 100644 index 00000000..01db5d5a --- /dev/null +++ b/crates/raw/tests/dng.rs @@ -0,0 +1,291 @@ +//! Synthetic DNG files: every storage layout decodes to the exact sensor +//! samples, and development reproduces known colours. + +use photocraft_raw::testgen::{DngSpec, DngStorage, mosaic, scene}; +use photocraft_raw::*; + +fn sensor(bytes: &[u8]) -> Sensor { + decode(bytes, &Limits::default()).unwrap() +} + +/// Deterministic sensor-like samples below `2^bits`. +fn noise(n: usize, bits: u16) -> Vec { + let mut s = 0x9E37_79B9_7F4A_7C15u64; + let max = (1u32 << bits) - 1; + (0..n) + .map(|i| { + s = s.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407); + ((i as u32 * 13 + (s >> 40) as u32 % 200) % (max + 1)) as u16 + }) + .collect() +} + +#[test] +fn identifies_dng() { + let b = DngSpec::cfa(4, 4, vec![0; 16]).build(); + assert_eq!(identify(&b), Some(RawFormat::Dng)); + assert!(is_raw(&b)); + assert!(!is_raw(b"II*\0\x08\0\0\0\0\0")); + assert!(!is_raw(b"\x89PNG\r\n\x1a\n")); +} + +#[test] +fn storage_layouts_decode_exactly() { + let (w, h) = (37, 23); // odd sizes exercise partial strips and edge tiles + let cases: Vec<(&str, u16, DngStorage, bool)> = vec![ + ("16-bit strips", 16, DngStorage::Strips { rows: 5 }, false), + ("16-bit big-endian", 16, DngStorage::Strips { rows: 23 }, true), + ("12-bit packed", 12, DngStorage::Strips { rows: 4 }, false), + ("10-bit packed big-endian", 10, DngStorage::Strips { rows: 7 }, true), + ("8-bit", 8, DngStorage::Strips { rows: 9 }, false), + ("LJ92 tiles", 14, DngStorage::Lj92Tiles { width: 16, height: 8 }, false), + ("LJ92 tiles big-endian", 16, DngStorage::Lj92Tiles { width: 32, height: 32 }, true), + ("LJ92 strips", 12, DngStorage::Lj92Strips { rows: 6 }, false), + ]; + for (name, bits, storage, be) in cases { + let data = noise(w * h, bits); + let mut spec = DngSpec::cfa(w, h, data.clone()); + spec.bits = bits; + spec.storage = storage; + spec.big_endian = be; + spec.white = (1 << bits) - 1; + let s = sensor(&spec.build()); + assert_eq!((s.width, s.height, s.samples), (w, h, 1), "{name}"); + assert_eq!(s.data, data, "{name}"); + assert_eq!(s.format, RawFormat::Dng); + } +} + +#[test] +fn linear_raw_decodes_and_develops() { + let (w, h) = (20, 10); + let data = noise(w * h * 3, 16); + let mut spec = DngSpec::cfa(w, h, data.clone()); + spec.samples = 3; + spec.storage = DngStorage::Lj92Tiles { width: 16, height: 16 }; + let s = sensor(&spec.build()); + assert_eq!((s.samples, s.cfa.is_none()), (3, true)); + assert_eq!(s.data, data); + let d = develop_sensor(&s, &DevelopOptions::default()).unwrap(); + assert_eq!((d.width, d.height, d.rgb.len()), (20, 10, 600)); + assert!(d.info.cfa.is_none()); +} + +#[test] +fn metadata_is_read() { + let (w, h) = (40, 30); + let mut spec = DngSpec::cfa(w, h, vec![1000; w * h]); + spec.bits = 14; + spec.white = 15000; + spec.black = vec![100, 101, 102, 103]; + spec.black_repeat = (2, 2); + spec.active_area = Some([2, 4, 28, 36]); // top, left, bottom, right + spec.default_crop = Some(([2, 1], [20, 16])); + spec.cfa = [1, 0, 2, 1]; // GRBG at the active-area origin + spec.as_shot_neutral = Some([0.5, 1.0, 0.8]); + spec.baseline_exposure = Some(0.5); + spec.orientation = 6; + spec.linearization = Some((0..4096u32).map(|v| (v * 3).min(65535) as u16).collect()); + let s = sensor(&spec.build()); + assert_eq!(s.white, [15000.0; 3]); + assert_eq!(s.black.values, vec![100.0, 101.0, 102.0, 103.0]); + assert_eq!(s.active, Rect::new(4, 2, 32, 26)); + assert_eq!(s.crop, Rect::new(6, 3, 20, 16)); + // Pattern anchored at the active area origin (4, 2): GRBG there. + assert_eq!(s.cfa.as_ref().unwrap().phase(4, 2), [1, 0, 2, 1]); + assert_eq!(s.color.as_shot_neutral, Some([0.5, 1.0, 0.8])); + assert!((s.baseline_exposure - 0.5).abs() < 1e-6); + assert_eq!(s.orientation, 6); + assert_eq!(s.linearization.as_ref().map(Vec::len), Some(4096)); + let d = develop_sensor(&s, &DevelopOptions::default()).unwrap(); + // Orientation 6 swaps the cropped 20×16 to 16×20. + assert_eq!((d.width, d.height), (16, 20)); + assert_eq!(d.info.wb_multipliers, [2.0, 1.0, 1.25]); + let d = develop_sensor(&s, &DevelopOptions { orient: false, ..Default::default() }).unwrap(); + assert_eq!((d.width, d.height), (20, 16)); +} + +/// Linear sRGB (D65) → XYZ, from the sRGB primaries (IEC 61966-2-1). +const SRGB_TO_XYZ: [[f64; 3]; 3] = [[0.4124564, 0.3575761, 0.1804375], [0.2126729, 0.7151522, 0.0721750], [0.0193339, 0.1191920, 0.9503041]]; + +fn invert(m: [[f64; 3]; 3]) -> [[f64; 3]; 3] { + let [[a, b, c], [d, e, f], [g, h, i]] = m; + let det = a * (e * i - f * h) - b * (d * i - f * g) + c * (d * h - e * g); + [ + [(e * i - f * h) / det, (c * h - b * i) / det, (b * f - c * e) / det], + [(f * g - d * i) / det, (a * i - c * g) / det, (c * d - a * f) / det], + [(d * h - e * g) / det, (b * g - a * h) / det, (a * e - b * d) / det], + ] +} + +fn ap(m: &[[f64; 3]; 3], v: [f64; 3]) -> [f64; 3] { + [0, 1, 2].map(|i| m[i][0] * v[0] + m[i][1] * v[1] + m[i][2] * v[2]) +} + +/// The DNG model computed independently: camera XYZ of `camera` under an +/// illuminant whose camera neutral is `neutral`, Bradford-adapted from that +/// white (normalized to Y = 1) to D50, then XYZ D50 → linear ProPhoto. +fn reference_prophoto(camera: [f64; 3], neutral: [f64; 3]) -> [f64; 3] { + const BRADFORD: [[f64; 3]; 3] = [[0.8951, 0.2664, -0.1614], [-0.7502, 1.7135, 0.0367], [0.0389, -0.0685, 1.0296]]; + const D50: [f64; 3] = [0.9642, 1.0, 0.8249]; + const XYZ_D50_TO_PROPHOTO: [[f64; 3]; 3] = [[1.3459433, -0.2556075, -0.0511118], [-0.5445989, 1.5081673, 0.0205351], [0.0, 0.0, 1.2118128]]; + let white = ap(&SRGB_TO_XYZ, neutral); + let white = white.map(|v| v / white[1]); + let xyz = ap(&SRGB_TO_XYZ, camera).map(|v| v / ap(&SRGB_TO_XYZ, neutral)[1]); + let (s, d) = (ap(&BRADFORD, white), ap(&BRADFORD, D50)); + let cone = ap(&BRADFORD, xyz); + let adapted = ap(&invert(BRADFORD), [cone[0] * d[0] / s[0], cone[1] * d[1] / s[1], cone[2] * d[2] / s[2]]); + ap(&XYZ_D50_TO_PROPHOTO, adapted) +} + +#[test] +fn colour_matrix_development_matches_reference() { + // A camera whose native RGB is linear sRGB: ColorMatrix = XYZ → sRGB at D65. + // The scene is captured with per-channel gains (the "white balance" the camera + // must undo), recorded as AsShotNeutral. + let (w, h) = (48, 32); + let gains = [0.5f32, 1.0, 0.7]; + let patches: [[f32; 3]; 4] = [[0.18, 0.18, 0.18], [0.40, 0.10, 0.05], [0.05, 0.30, 0.10], [0.10, 0.12, 0.45]]; + let rgb: Vec<[f32; 3]> = (0..w * h) + .map(|i| { + let p = patches[((i % w) / (w / 2)) + 2 * ((i / w) / (h / 2))]; + [p[0] * gains[0], p[1] * gains[1], p[2] * gains[2]] + }) + .collect(); + let mut spec = DngSpec::cfa(w, h, mosaic(&rgb, w, [0, 1, 1, 2], 256, 16383)); + spec.bits = 14; + spec.black = vec![256]; + spec.white = 16383; + let cm = invert(SRGB_TO_XYZ); + spec.color_matrix1 = Some((21, [cm[0][0], cm[0][1], cm[0][2], cm[1][0], cm[1][1], cm[1][2], cm[2][0], cm[2][1], cm[2][2]])); + spec.as_shot_neutral = Some(gains.map(f64::from)); + for storage in [DngStorage::Strips { rows: 32 }, DngStorage::Lj92Tiles { width: 16, height: 16 }] { + spec.storage = storage; + for method in [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] { + let d = develop(&spec.build(), &DevelopOptions { demosaic: method, ..Default::default() }).unwrap(); + for (pi, p) in patches.iter().enumerate() { + // Centre of each patch. + let (cx, cy) = ((pi % 2) * (w / 2) + w / 4, (pi / 2) * (h / 2) + h / 4); + let o = &d.rgb[(cy * w + cx) * 3..(cy * w + cx) * 3 + 3]; + let g = gains.map(f64::from); + let want = reference_prophoto([0, 1, 2].map(|c| f64::from(p[c]) * g[c]), g); + for c in 0..3 { + let got = (f64::from(o[c]) / 65535.0).powf(1.8); + assert!((got - want[c]).abs() < 0.006, "{storage:?} {method:?} patch {pi} c {c}: got {got:.4} want {:.4}", want[c]); + } + } + } + } +} + +#[test] +fn forward_matrix_path_maps_neutral_to_white() { + let (w, h) = (16, 16); + let rgb = vec![[0.3f32, 0.6, 0.45]; w * h]; // a grey under a green-ish light + let mut spec = DngSpec::cfa(w, h, mosaic(&rgb, w, [0, 1, 1, 2], 0, 65535)); + let cm = invert(SRGB_TO_XYZ); + let flat = |m: [[f64; 3]; 3]| [m[0][0], m[0][1], m[0][2], m[1][0], m[1][1], m[1][2], m[2][0], m[2][1], m[2][2]]; + spec.color_matrix1 = Some((17, flat(cm))); + spec.color_matrix2 = Some((21, flat(cm))); + // ForwardMatrix: balanced camera → XYZ D50 (sRGB adapted to D50; white → D50). + let fm = [[0.4360747, 0.3850649, 0.1430804], [0.2225045, 0.7168786, 0.0606169], [0.0139322, 0.0971045, 0.7141733]]; + spec.forward_matrix1 = Some(flat(fm)); + spec.forward_matrix2 = Some(flat(fm)); + spec.as_shot_neutral = Some([0.5, 1.0, 0.75]); + let d = develop(&spec.build(), &DevelopOptions::default()).unwrap(); + let o = &d.rgb[(8 * w + 8) * 3..(8 * w + 8) * 3 + 3]; + // Balanced grey 0.6 → ProPhoto (0.6, 0.6, 0.6). + for (c, v) in o.iter().enumerate() { + let got = (f64::from(*v) / 65535.0).powf(1.8); + assert!((got - 0.6).abs() < 0.004, "c {c}: {got}"); + } +} + +#[test] +fn highlights_clip_to_white_and_exposure_applies() { + let (w, h) = (8, 8); + let mut spec = DngSpec::cfa(w, h, vec![65535; w * h]); + spec.as_shot_neutral = Some([0.5, 1.0, 0.6]); + let d = develop(&spec.build(), &DevelopOptions::default()).unwrap(); + assert!(d.rgb.iter().all(|&v| v >= 65500), "clipped sensor must render white"); + let mut spec = DngSpec::cfa(w, h, vec![16384; w * h]); // 0.25 linear grey + spec.as_shot_neutral = Some([1.0, 1.0, 1.0]); + let b = spec.build(); + let base = develop(&b, &DevelopOptions::default()).unwrap().rgb[27]; + let up = develop(&b, &DevelopOptions { exposure: 1.0, ..Default::default() }).unwrap().rgb[27]; + let lin = |v: u16| (f64::from(v) / 65535.0).powf(1.8); + assert!((lin(up) / lin(base) - 2.0).abs() < 0.02, "{} vs {}", lin(up), lin(base)); +} + +#[test] +fn gain_map_opcode_is_applied() { + // Uniform 0.25 grey; a gain map that doubles everything over the whole image. + let (w, h) = (16, 12); + let mut spec = DngSpec::cfa(w, h, vec![16384; w * h]); + spec.as_shot_neutral = Some([1.0; 3]); + let lin = |v: u16| (f64::from(v) / 65535.0).powf(1.8); + let base = develop(&spec.build(), &DevelopOptions::default()).unwrap().rgb[(6 * w + 8) * 3 + 1]; + spec.opcode_list2 = Some(photocraft_raw::testgen::gain_map_opcode_list([0, 0, h as u32, w as u32], 1, [[2.0, 2.0], [2.0, 2.0]])); + let d = develop(&spec.build(), &DevelopOptions::default()).unwrap(); + assert!(!d.warnings.iter().any(|w| w.contains("Opcode")), "{:?}", d.warnings); + let gained = d.rgb[(6 * w + 8) * 3 + 1]; + assert!((lin(gained) / lin(base) - 2.0).abs() < 0.02, "{} vs {}", lin(gained), lin(base)); +} + +#[test] +fn white_balance_options() { + let (w, h) = (16, 16); + let rgb = vec![[0.2f32, 0.4, 0.3]; w * h]; + let spec = DngSpec::cfa(w, h, mosaic(&rgb, w, [0, 1, 1, 2], 0, 65535)); + let b = spec.build(); + // No AsShotNeutral: grey world → neutral output. + let d = develop(&b, &DevelopOptions::default()).unwrap(); + assert!(d.warnings.iter().any(|w| w.contains("grey world"))); + let o = &d.rgb[(8 * w + 8) * 3..(8 * w + 8) * 3 + 3]; + assert!((i32::from(o[0]) - i32::from(o[1])).abs() < 300 && (i32::from(o[2]) - i32::from(o[1])).abs() < 300, "{o:?}"); + let d = develop(&b, &DevelopOptions { white_balance: WhiteBalance::Multipliers([2.0, 1.0, 4.0 / 3.0]), ..Default::default() }).unwrap(); + assert_eq!(d.info.wb_multipliers.map(|v| (v * 1000.0).round()), [2000.0, 1000.0, 1333.0]); + assert!(develop(&b, &DevelopOptions { white_balance: WhiteBalance::Multipliers([0.0, 1.0, 1.0]), ..Default::default() }).is_err()); +} + +#[test] +fn smooth_scene_develops_with_every_method() { + let (w, h) = (64, 48); + let rgb = scene(w, h); + let spec = DngSpec::cfa(w, h, mosaic(&rgb, w, [2, 1, 1, 0], 0, 65535)); + let mut spec = spec; + spec.cfa = [2, 1, 1, 0]; + spec.as_shot_neutral = Some([1.0; 3]); + let b = spec.build(); + let outs: Vec> = [Demosaic::Bilinear, Demosaic::Mhc, Demosaic::Ahd] + .iter() + .map(|&m| develop(&b, &DevelopOptions { demosaic: m, ..Default::default() }).unwrap().rgb) + .collect(); + // Interior pixels agree closely across methods on a smooth scene. + for y in 4..h - 4 { + for x in 4..w - 4 { + for c in 0..3 { + let i = (y * w + x) * 3 + c; + let (a, b, d) = (i32::from(outs[0][i]), i32::from(outs[1][i]), i32::from(outs[2][i])); + assert!((a - b).abs() < 800 && (a - d).abs() < 800, "({x},{y}) c{c}: {a} {b} {d}"); + } + } + } +} + +#[test] +fn unsupported_variants_say_so() { + let mut spec = DngSpec::cfa(8, 8, vec![0; 64]); + spec.cfa = [0, 1, 2, 1]; // not Bayer + assert!(matches!(decode(&spec.build(), &Limits::default()), Err(RawError::Unsupported(_)))); +} + +#[test] +fn limits_are_enforced_before_allocating() { + let spec = DngSpec::cfa(64, 64, vec![0; 64 * 64]); + let b = spec.build(); + let tiny = Limits { max_width: 32, ..Limits::default() }; + assert!(matches!(decode(&b, &tiny), Err(RawError::LimitExceeded(_)))); + let tiny = Limits { max_alloc: 1000, ..Limits::default() }; + assert!(matches!(decode(&b, &tiny), Err(RawError::LimitExceeded(_)))); +} diff --git a/crates/raw/tests/hostile.rs b/crates/raw/tests/hostile.rs new file mode 100644 index 00000000..e78ee337 --- /dev/null +++ b/crates/raw/tests/hostile.rs @@ -0,0 +1,115 @@ +//! Truncated and corrupted raw files must fail cleanly, never panic, and +//! never allocate beyond the limits. + +use photocraft_raw::testgen::{Cr2Spec, DngSpec, DngStorage, mosaic, scene, tiff_ep}; +use photocraft_raw::*; + +fn samples() -> Vec> { + let (w, h) = (24, 12); + let cfa = mosaic(&scene(w, h), w, [0, 1, 1, 2], 64, 4000); + let mut out = Vec::new(); + let mut d = DngSpec::cfa(w, h, cfa.clone()); + d.bits = 12; + d.white = 4000; + d.black = vec![64]; + d.as_shot_neutral = Some([0.5, 1.0, 0.7]); + d.color_matrix1 = Some((21, [3.2, -1.5, -0.5, -1.0, 1.9, 0.04, 0.06, -0.2, 1.06])); + out.push(d.build()); + d.storage = DngStorage::Lj92Tiles { width: 16, height: 8 }; + d.active_area = Some([1, 2, 11, 22]); + d.default_crop = Some(([1, 1], [16, 8])); + out.push(d.build()); + d.storage = DngStorage::Strips { rows: 5 }; + d.big_endian = true; + out.push(d.build()); + let mut data = vec![600u16; 40 * 10]; + for (i, v) in data.iter_mut().enumerate() { + *v = 600 + (i as u16 * 7) % 3000; + } + out.push( + Cr2Spec { + width: 40, + height: 10, + data, + precision: 14, + components: 2, + slices: vec![16, 16, 8], + borders: Some([8, 2, 39, 9]), + wb_rggb: Some([2000, 1024, 1024, 1500]), + orientation: 8, + } + .build(), + ); + out.push(tiff_ep("NIKON", w, h, &cfa, [0, 1, 1, 2], 12, vec![])); + out +} + +fn exercise(b: &[u8]) { + let _ = identify(b); + let _ = embedded_preview(b); + let opts = DevelopOptions { limits: Limits { max_alloc: 64 << 20, ..Limits::default() }, ..Default::default() }; + if let Ok(s) = decode(b, &opts.limits) { + for m in [Demosaic::Bilinear, Demosaic::Ahd] { + let _ = develop_sensor(&s, &DevelopOptions { demosaic: m, ..opts.clone() }); + } + } +} + +#[test] +fn originals_develop() { + for (i, b) in samples().iter().enumerate() { + let d = develop(b, &DevelopOptions::default()); + assert!(d.is_ok(), "sample {i}: {:?}", d.err()); + } +} + +#[test] +fn every_truncation_fails_cleanly() { + for b in samples() { + for n in (0..b.len()).step_by(3) { + exercise(&b[..n]); + } + } +} + +#[test] +fn random_corruption_never_panics() { + let mut s = 0xDEAD_BEEF_u64; + let mut next = || { + s = s.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407); + s >> 33 + }; + for b in samples() { + for _ in 0..400 { + let mut c = b.clone(); + for _ in 0..1 + next() % 6 { + let at = next() as usize % c.len(); + c[at] = match next() % 4 { + 0 => 0xFF, + 1 => 0, + _ => next() as u8, + }; + } + exercise(&c); + } + } +} + +#[test] +fn absurd_dimensions_are_rejected_without_allocating() { + let mut d = DngSpec::cfa(4, 4, vec![0; 16]); + d.storage = DngStorage::Strips { rows: 4 }; + let mut b = d.build(); + // Patch the raw IFD's ImageWidth / ImageLength (LONG, inline) to 200000 × 200000. + let t = 256u16.to_le_bytes(); + let mut patched = 0; + for i in 0..b.len().saturating_sub(12) { + if b[i..i + 2] == t && b[i + 2..i + 4] == 4u16.to_le_bytes() && b[i + 8..i + 12] == 4u32.to_le_bytes() { + b[i + 8..i + 12].copy_from_slice(&200_000u32.to_le_bytes()); + b[i + 20..i + 24].copy_from_slice(&200_000u32.to_le_bytes()); + patched += 1; + } + } + assert!(patched > 0); + assert!(matches!(decode(&b, &Limits::default()), Err(RawError::LimitExceeded(_)))); +} diff --git a/crates/raw/tests/vendors.rs b/crates/raw/tests/vendors.rs new file mode 100644 index 00000000..561958a4 --- /dev/null +++ b/crates/raw/tests/vendors.rs @@ -0,0 +1,132 @@ +//! Synthetic CR2 and TIFF/EP (NEF / ARW-like) files. + +use photocraft_raw::testgen::{Cr2Spec, Val, mosaic, scene, tiff_ep}; +use photocraft_raw::*; + +fn cr2_spec(w: usize, h: usize, comps: usize, slices: Vec) -> Cr2Spec { + // Left / top borders of 16 / 4 masked pixels at black level 512. + let (l, t) = (16, 4); + let rgb = scene(w - l, h - t); + let img = mosaic(&rgb, w - l, [0, 1, 1, 2], 512, 13000); + let mut data = vec![512u16; w * h]; + for y in t..h { + for x in l..w { + data[y * w + x] = img[(y - t) * (w - l) + (x - l)]; + } + } + Cr2Spec { + width: w, + height: h, + data, + precision: 14, + components: comps, + slices, + borders: Some([l as u16, t as u16, w as u16 - 1, h as u16 - 1]), + wb_rggb: Some([2048, 1024, 1024, 1536]), + orientation: 1, + } +} + +#[test] +fn cr2_slices_decode_exactly() { + for (comps, slices) in [(2, vec![32, 32, 24]), (4, vec![40, 48]), (2, vec![]), (4, vec![88])] { + let spec = cr2_spec(88, 30, comps, slices.clone()); + let b = spec.build(); + assert_eq!(identify(&b), Some(RawFormat::Cr2)); + let s = decode(&b, &Limits::default()).unwrap(); + assert_eq!((s.width, s.height), (88, 30), "{comps} {slices:?}"); + assert_eq!(s.data, spec.data, "{comps} components, slices {slices:?}"); + assert_eq!(s.active, Rect::new(16, 4, 72, 26)); + assert_eq!(s.camera_wb, Some([2.0, 1.0, 1.5])); + assert!(s.black.values.iter().all(|&v| (v - 512.0).abs() < 1e-3), "{:?}", s.black.values); + assert_eq!(s.cfa.as_ref().unwrap().phase(16, 4), [0, 1, 1, 2]); + } +} + +#[test] +fn cr2_develops() { + let b = cr2_spec(88, 30, 2, vec![32, 32, 24]).build(); + let d = develop(&b, &DevelopOptions::default()).unwrap(); + assert_eq!((d.width, d.height), (72, 26)); + assert_eq!(d.info.format, RawFormat::Cr2); + assert_eq!(d.info.wb_multipliers, [2.0, 1.0, 1.5]); + assert!(d.warnings.iter().any(|w| w.contains("sRGB primaries"))); +} + +#[test] +fn cr2_srgb_raw_is_unsupported() { + // A subsampled (sRAW-style) frame: patch the first component's sampling factors. + let mut b = cr2_spec(32, 8, 2, vec![]).build(); + let sof = b.windows(2).position(|w| w == [0xFF, 0xC3]).unwrap(); + b[sof + 11] = 0x21; + match decode(&b, &Limits::default()) { + Err(RawError::Unsupported(m)) => assert!(m.contains("sRAW"), "{m}"), + other => panic!("expected unsupported sRAW, got {other:?}"), + } +} + +#[test] +fn nef_like_with_maker_tags() { + let (w, h) = (24, 16); + let data = mosaic(&scene(w, h), w, [2, 1, 1, 0], 0, 4095); + let b = tiff_ep("NIKON CORPORATION", w, h, &data, [2, 1, 1, 0], 12, vec![]); + assert_eq!(identify(&b), Some(RawFormat::Nef)); + let s = decode(&b, &Limits::default()).unwrap(); + assert_eq!(s.data, data); + assert_eq!(s.cfa.as_ref().unwrap().phase(0, 0), [2, 1, 1, 0]); + assert!(s.warnings.iter().any(|w| w.contains("black level"))); + assert!(develop_sensor(&s, &DevelopOptions::default()).is_ok()); +} + +#[test] +fn arw_like_reads_sony_levels_and_crop() { + let (w, h) = (24, 16); + let data = mosaic(&scene(w, h), w, [0, 1, 1, 2], 512, 16383); + let extra = vec![ + (0x7310, Val::Short(vec![512, 512, 512, 512])), + (0x7313, Val::Short(vec![2048, 1024, 1024, 1536])), + (50717, Val::Short(vec![16383])), + (50719, Val::Long(vec![2, 2])), + (50720, Val::Long(vec![20, 12])), + ]; + let b = tiff_ep("SONY", w, h, &data, [0, 1, 1, 2], 14, extra); + assert_eq!(identify(&b), Some(RawFormat::Arw)); + let s = decode(&b, &Limits::default()).unwrap(); + assert_eq!(s.black.values, vec![512.0; 4]); + assert_eq!(s.white, [16383.0; 3]); + assert_eq!(s.camera_wb, Some([2.0, 1.0, 1.5])); + assert_eq!(s.crop, Rect::new(2, 2, 20, 12)); + assert!(s.warnings.is_empty(), "{:?}", s.warnings); + let d = develop_sensor(&s, &DevelopOptions::default()).unwrap(); + assert_eq!((d.width, d.height), (20, 12)); +} + +#[test] +fn other_raw_containers_are_recognised_but_unsupported() { + let mut cr3 = vec![0, 0, 0, 24]; + cr3.extend_from_slice(b"ftypcrx "); + cr3.extend_from_slice(&[0; 12]); + assert_eq!(identify(&cr3), Some(RawFormat::Cr3)); + assert!(matches!(decode(&cr3, &Limits::default()), Err(RawError::Unsupported(_)))); + let mut raf = b"FUJIFILMCCD-RAW 0201FF383501".to_vec(); + raf.resize(200, 0); + assert_eq!(identify(&raf), Some(RawFormat::Raf)); + assert!(matches!(decode(&raf, &Limits::default()), Err(RawError::Unsupported(_)))); + assert_eq!(decode(b"not a raw file", &Limits::default()).unwrap_err(), RawError::NotRaw); +} + +#[test] +fn embedded_preview_is_found() { + // A RAF header pointing at a minimal baseline JPEG (SOI, SOF0, EOI). + let jpeg = [0xFF, 0xD8, 0xFF, 0xC0, 0x00, 0x0B, 8, 0x00, 0x20, 0x00, 0x30, 1, 1, 0x11, 0, 0xFF, 0xD9]; + let mut raf = b"FUJIFILMCCD-RAW 0201FF383501".to_vec(); + raf.resize(100, 0); + raf[84..88].copy_from_slice(&100u32.to_be_bytes()); + raf[88..92].copy_from_slice(&(jpeg.len() as u32).to_be_bytes()); + raf.extend_from_slice(&jpeg); + let p = embedded_preview(&raf).unwrap(); + assert_eq!((p.width, p.height), (0x30, 0x20)); + assert_eq!(p.jpeg, &jpeg); + // The lossless-JPEG raw data of a CR2 is not mistaken for a preview. + assert!(embedded_preview(&cr2_spec(32, 8, 2, vec![]).build()).is_none()); +} diff --git a/docs/architecture.md b/docs/architecture.md index da7eb37d..26f6d536 100644 --- a/docs/architecture.md +++ b/docs/architecture.md @@ -98,7 +98,7 @@ photocraft/ │ ├─ psd/ photocraft-psd PSD/PSB read + write; its OWN format-level model; depends on nothing in this workspace │ ├─ adobe-assets/ photocraft-adobe-assets .abr .asl .aco/.ase .grd .pat .csh .atn .cube/.3dl, ACR .xmp presets (standalone, like psd) │ ├─ codecs/ photocraft-codecs png/jpeg/tiff/webp/gif/bmp/avif/jxl (+heif via feature) decode/encode -│ ├─ raw/ photocraft-raw camera RAW decode (rawler) + develop pipeline + lens profiles +│ ├─ raw/ photocraft-raw clean-room camera RAW decode (DNG, CR2, TIFF/EP) + develop pipeline (standalone, like psd) │ ├─ format/ photocraft-format native document format (.pcraft bundle): manifest + content-addressed tiles │ ├─ io/ photocraft-io import/export orchestration; doc ⇄ PSD mapping; PDF/SVG import (features) │ │ ── intelligence ── @@ -122,7 +122,7 @@ photocraft/ └─ plan/ this directory ``` -**What exists today.** This layout is the target design. Built so far: `geom`, `cms`, `color`, `raster`, `psd`, `codecs`, `doc`, `ops`, `algo`, `paint`, `text`, `vector`, `compose`, `gpu`, `format`, `io`, `engine`, `ui-egui`, `automation`, `testkit`, and the three apps. Not yet split out: `viewport` and `tools` live inside `ui-egui` and `engine`; `platform` services are function hooks injected by each app (`ui_egui::Services`); `adobe-assets`, `raw` and `ml` are not started. +**What exists today.** This layout is the target design. Built so far: `geom`, `cms`, `color`, `raster`, `psd`, `codecs`, `doc`, `ops`, `algo`, `paint`, `text`, `vector`, `compose`, `gpu`, `format`, `raw`, `io`, `engine`, `ui-egui`, `automation`, `testkit`, and the three apps. Not yet split out: `viewport` and `tools` live inside `ui-egui` and `engine`; `platform` services are function hooks injected by each app (`ui_egui::Services`); `adobe-assets` and `ml` are not started. **Crate granularity:** start with the crates above. Split `algo` into `-adjust`, `-filters`, `-select`, `-inpaint` and `-warp` once any module passes about 10k lines, or once compile times hurt. Its internal module boundaries should already follow those lines. @@ -137,10 +137,10 @@ photocraft/ L6 ui-egui · automation · platform L5 engine L4 tools · viewport · io · ml - L3 compose · gpu · format · raw · codecs + L3 compose · gpu · format L2 ops · algo · paint · text · vector L1 doc - L0 geom · color · raster psd, adobe-assets (standalone, no workspace deps) + L0 geom · cms · color · raster psd, codecs, raw, adobe-assets (standalone, no workspace deps) ``` **Rules** (checked by `cargo xtask layers` in CI, which parses `cargo metadata`): diff --git a/xtask/src/layers.rs b/xtask/src/layers.rs index d39b7a4b..6dd7e0be 100644 --- a/xtask/src/layers.rs +++ b/xtask/src/layers.rs @@ -39,6 +39,7 @@ pub const TABLE: &[(&str, Class)] = &[ ("raster", Class::Layer(0)), ("psd", Class::Standalone), ("codecs", Class::Standalone), + ("raw", Class::Standalone), ("adobe-assets", Class::Standalone), ("doc", Class::Layer(1)), ("ops", Class::Layer(2)), @@ -49,7 +50,6 @@ pub const TABLE: &[(&str, Class)] = &[ ("compose", Class::Layer(3)), ("gpu", Class::Layer(3)), ("format", Class::Layer(3)), - ("raw", Class::Layer(3)), ("io", Class::Layer(4)), ("tools", Class::Layer(4)), ("viewport", Class::Layer(4)),