diff --git a/crates/j2k-native/fixtures/subsampled-origin.j2c b/crates/j2k-native/fixtures/subsampled-origin.j2c
new file mode 100644
index 00000000..08682df7
Binary files /dev/null and b/crates/j2k-native/fixtures/subsampled-origin.j2c differ
diff --git a/crates/j2k-native/src/j2c/decode/store.rs b/crates/j2k-native/src/j2c/decode/store.rs
index af6c9035..61b29cfc 100644
--- a/crates/j2k-native/src/j2c/decode/store.rs
+++ b/crates/j2k-native/src/j2c/decode/store.rs
@@ -189,13 +189,15 @@ pub(super) fn store<'a>(
 
         let (x_offset, y_offset) = output_grid_offset(&header.size_data);
 
-        // Otherwise, copy sample by sample.
-        for y in resolution_tile.rect.y0..resolution_tile.rect.y1 {
-            let relative_y = (y - component_tile.rect.y0) as usize;
+        // Coefficients use the reconstructed resolution grid, which may be
+        // reduced or clipped relative to the full component tile.
+        let sample_rect = resolution_tile.rect.intersect(idwt_output.rect);
+        for y in sample_rect.y0..sample_rect.y1 {
+            let relative_y = (y - idwt_output.rect.y0) as usize;
             let reference_grid_y = (u32::from(scale_y) * y) / y_shrink_factor;
 
-            for x in resolution_tile.rect.x0..resolution_tile.rect.x1 {
-                let relative_x = (x - component_tile.rect.x0) as usize;
+            for x in sample_rect.x0..sample_rect.x1 {
+                let relative_x = (x - idwt_output.rect.x0) as usize;
                 let reference_grid_x = (u32::from(scale_x) * x) / x_shrink_factor;
 
                 let sample = idwt_output.coefficients
@@ -327,12 +329,13 @@ fn store_i64<'a>(
 
         let (x_offset, y_offset) = output_grid_offset(&header.size_data);
 
-        for y in resolution_tile.rect.y0..resolution_tile.rect.y1 {
-            let relative_y = (y - component_tile.rect.y0) as usize;
+        let sample_rect = resolution_tile.rect.intersect(idwt_output.rect);
+        for y in sample_rect.y0..sample_rect.y1 {
+            let relative_y = (y - idwt_output.rect.y0) as usize;
             let reference_grid_y = (u32::from(scale_y) * y) / y_shrink_factor;
 
-            for x in resolution_tile.rect.x0..resolution_tile.rect.x1 {
-                let relative_x = (x - component_tile.rect.x0) as usize;
+            for x in sample_rect.x0..sample_rect.x1 {
+                let relative_x = (x - idwt_output.rect.x0) as usize;
                 let reference_grid_x = (u32::from(scale_x) * x) / x_shrink_factor;
 
                 let sample = idwt_output.coefficients_i64
@@ -516,12 +519,13 @@ fn store_region<'a>(
         round_irreversible_output,
     );
 
-    for y in resolution_tile.rect.y0..resolution_tile.rect.y1 {
-        let relative_y = (y - component_tile.rect.y0) as usize;
+    let sample_rect = resolution_tile.rect.intersect(idwt_output.rect);
+    for y in sample_rect.y0..sample_rect.y1 {
+        let relative_y = (y - idwt_output.rect.y0) as usize;
         let reference_grid_y = (u32::from(scale_y) * y) / y_shrink_factor;
 
-        for x in resolution_tile.rect.x0..resolution_tile.rect.x1 {
-            let relative_x = (x - component_tile.rect.x0) as usize;
+        for x in sample_rect.x0..sample_rect.x1 {
+            let relative_x = (x - idwt_output.rect.x0) as usize;
             let reference_grid_x = (u32::from(scale_x) * x) / x_shrink_factor;
 
             let sample = idwt_output.coefficients
diff --git a/crates/j2k-native/tests/component_planes.rs b/crates/j2k-native/tests/component_planes.rs
index 6656972b..ff368981 100644
--- a/crates/j2k-native/tests/component_planes.rs
+++ b/crates/j2k-native/tests/component_planes.rs
@@ -866,3 +866,32 @@ fn classic_reversible_i64_decode_native_region_crops_29_bit_bytes() {
     assert_eq!(components.planes()[0].bit_depth(), 29);
     assert_eq!(components.planes()[0].data(), expected);
 }
+
+#[test]
+fn subsampled_nonzero_origin_decodes_full_and_reduced_regions() {
+    // Minimized wsi-rs open_jp2k_codestream_bytes fuzz regression: XOsiz=2,
+    // and component 1 has vertical sampling 47 with two decomposition levels.
+    let bytes = include_bytes!("../fixtures/subsampled-origin.j2c");
+    for reduction in 0..=2 {
+        let image = Image::new_with_reduction(bytes, &DecodeSettings::strict(), reduction)
+            .expect("parse subsampled nonzero-origin codestream");
+        let region = image
+            .decode_region((0, 0, image.width(), image.height()))
+            .expect("decode full image region without coordinate underflow");
+        let full = image.decode().expect("decode full subsampled image");
+        assert_eq!(region.data, full, "reduction {reduction}");
+
+        let width = image.width() as usize;
+        let height = image.height() as usize;
+        let cropped = image
+            .decode_region((1, 1, image.width() - 2, image.height() - 2))
+            .expect("decode clipped subsampled region");
+        let expected: Vec<u8> = full
+            .chunks_exact(width * 3)
+            .skip(1)
+            .take(height - 2)
+            .flat_map(|row| row[3..(width - 1) * 3].iter().copied())
+            .collect();
+        assert_eq!(cropped.data, expected, "cropped reduction {reduction}");
+    }
+}