diff --git a/benchmark.py b/benchmark.py
index 06e5f82..f448128 100644
--- a/benchmark.py
+++ b/benchmark.py
@@ -20,7 +20,7 @@ def main():
     # cal_path为校准文件路径, 每一个设备需要创建/使用一个校准文件
 
     # 方式1: 从视频文件读取
-    sensor = orisys.Sensor(0, config_name="./config/prototype_of3.json", cal_path= "./config/dd01a_0.npy", isstitch=True, cuda=False, verbose=True)
+    sensor = orisys.Sensor(2, config_name="./config/ddjx01_q6a.json", cal_path= "./config/ddjx01.npy", isstitch=True, cuda=False, verbose=True, backend='cpu')
 
     # 方式2: 从摄像头读取（取消注释使用） 
     # sensor = orisys.Sensor(0, config_name="./config/prototype_of3.json", cal_path= "./config/cal.npy", isstitch=True, cuda=True, verbose=True)
@@ -52,7 +52,7 @@ def main():
 
         # 步骤 2.3: 计算触碰点 
         start_time = time.perf_counter()
-        sensor.compute_contact()
+        # sensor.compute_contact()
         end_time = time.perf_counter()
         compute_contact_time = end_time - start_time
 
diff --git a/config/compound_eye.png b/config/compound_eye.png
index 117fba7..3b89059 100644
Binary files a/config/compound_eye.png and b/config/compound_eye.png differ
diff --git a/config/ddjx01.json b/config/ddjx01.json
index dfcf2f8..988ab1e 100644
--- a/config/ddjx01.json
+++ b/config/ddjx01.json
@@ -27,8 +27,9 @@
     "loop_closure_interval": 1,
     "algorithm": "dis",
     "dis_preset": "medium",
-    "dis_mode": 1,
-    "algorithm_list": ["raft", "dis"]
+    "dis_mode": 3,
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/config/ddjx01_480p.json b/config/ddjx01_480p.json
index ac85a9d..6083c90 100644
--- a/config/ddjx01_480p.json
+++ b/config/ddjx01_480p.json
@@ -11,8 +11,8 @@
   },
   "stitcher": {
     "min_area_frac": 0.01,
-    "max_area_frac": 0.06,
-    "max_aspect_ratio": 50,
+    "max_area_frac": 0.08,
+    "max_aspect_ratio": 100,
     "index": [1,2,3,6,7,8,11,12,13],
     "margin": 4,
     "nominal_ov_pos": 44,
@@ -28,7 +28,8 @@
     "algorithm": "dis",
     "dis_preset": "medium",
     "dis_mode": 2,
-    "algorithm_list": ["raft", "dis"]
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/config/ddjx01_q6a.json b/config/ddjx01_q6a.json
new file mode 100644
index 0000000..f267045
--- /dev/null
+++ b/config/ddjx01_q6a.json
@@ -0,0 +1,46 @@
+{
+    "name": "prototype",
+    "version": "1.0",
+    "camera": {
+      "resolution_raw": [1280, 720],
+      "fps": 30,
+      "units": 15,
+      "units_used": 9, 
+      "resolution_stitch": [263, 263],
+      "stitch_res_list": [365,263,123,63]
+    },
+    "stitcher": {
+      "min_area_frac": 0.01,
+      "max_area_frac": 0.04,
+      "max_aspect_ratio": 25,
+      "index": [1,2,3,6,7,8,11,12,13],
+      "margin": 10,
+      "nominal_ov_pos": 67,
+      "n": [3,3],
+      "row_tol": 0.6,
+      "box_size": [180, 180]
+    },
+    "tracker": {
+      "adaptive": true,
+      "lower_threshold": 12.0,
+      "upper_threshold": 18.0,
+      "loop_closure_interval": 1,
+      "algorithm": "dis",
+      "dis_preset": "medium",
+      "dis_mode": 3,
+      "algorithm_list": ["raft", "dis"],
+      "hhd_downsample": 3,
+      "backend": "cpu"
+    },
+    "contact": {
+      "min_contour_area": 3000,
+      "flow_threshold":1,
+      "div_threshold":0.02,
+      "expansion_threshold": 0.01,
+      "kernel_size": [5,5],
+      "denoise_sigma": 10
+    }
+  }
+  
+  
+  
\ No newline at end of file
diff --git a/config/prototype_of3.json b/config/prototype_of3.json
index dc1e174..8d94471 100644
--- a/config/prototype_of3.json
+++ b/config/prototype_of3.json
@@ -28,7 +28,8 @@
     "algorithm": "dis",
     "dis_preset": "medium",
     "dis_mode": 3,
-    "algorithm_list": ["raft", "dis"]
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/config/prototype_of3_480p.json b/config/prototype_of3_480p.json
index 52bb939..5c321a3 100644
--- a/config/prototype_of3_480p.json
+++ b/config/prototype_of3_480p.json
@@ -28,7 +28,8 @@
     "algorithm": "dis",
     "dis_preset": "medium",
     "dis_mode": 2,
-    "algorithm_list": ["raft", "dis"]
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/config/prototype_of3_720p.json b/config/prototype_of3_720p.json
index 3575674..185ccbf 100644
--- a/config/prototype_of3_720p.json
+++ b/config/prototype_of3_720p.json
@@ -18,7 +18,7 @@
     "nominal_ov_pos": 67,
     "n": [3,3],
     "row_tol": 0.6,
-    "box_size": [140, 140]
+    "box_size": [180, 180]
   },
   "tracker": {
     "adaptive": true,
@@ -28,7 +28,8 @@
     "algorithm": "dis",
     "dis_preset": "medium",
     "dis_mode": 2,
-    "algorithm_list": ["raft", "dis"]
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/config/prototype_raft.json b/config/prototype_raft.json
index 80883ee..480995b 100644
--- a/config/prototype_raft.json
+++ b/config/prototype_raft.json
@@ -27,7 +27,8 @@
     "algorithm": "raft",
     "dis_preset": "medium",
     "dis_mode": 3,
-    "algorithm_list": ["raft", "dis"]
+    "algorithm_list": ["raft", "dis"],
+    "hhd_downsample": 3
   },
   "contact": {
     "min_contour_area": 3000,
diff --git a/example2_two_sensor_pipeline.py b/example2_two_sensor_pipeline.py
index 09dc58f..60e5fd4 100644
--- a/example2_two_sensor_pipeline.py
+++ b/example2_two_sensor_pipeline.py
@@ -66,7 +66,7 @@ def sensor_worker(sensor_id, vid_src, config_name, cal_path, result_queue, stop_
                 sensor.compute_deformation()
 
                 # 步骤3: 计算触碰点
-                sensor.compute_contact()
+                # sensor.compute_contact()
 
                 # 步骤4: 读取信息
                 fps, flow, vnormal, img, contour, centroid, depth_map = sensor.read_info(
@@ -130,14 +130,14 @@ def main():
     sensor_configs = [
         {
             'sensor_id': 0,
-            'vid_src': 0,
-            'config_name': './config/ddjx01.json',
+            'vid_src': 2,
+            'config_name': './config/prototype_of3_720p.json',
             'cal_path': './config/ddjx01.npy'
         },  
         {
             'sensor_id': 1,
-            'vid_src': 1,  # 可以改为不同的视频源或摄像头索引
-            'config_name': './config/ddjx01.json',
+            'vid_src': 4,  # 可以改为不同的视频源或摄像头索引
+            'config_name': './config/prototype_of3_720p.json',
             'cal_path': './config/ddjx01.npy'
         },
     ]
diff --git a/src/orisys/HHD.py b/src/orisys/HHD.py
index 225f712..3c6f919 100644
--- a/src/orisys/HHD.py
+++ b/src/orisys/HHD.py
@@ -1626,7 +1626,8 @@ class PoissonSolver(object):
 
     def _get_opencl_batch_mul_kernel(self, queue):
         """Batched complex multiply: multiplies each batch plane by the same Green's spectrum.
-        Uses multiple elements per work item and explicit local size for better throughput."""
+        Uses multiple elements per work item, __constant for broadcast-read Green's spectrum,
+        and fma for better throughput/precision."""
         if self._opencl_batch_mul_kernel is not None:
             return self._opencl_batch_mul_kernel
 
@@ -1637,13 +1638,14 @@ class PoissonSolver(object):
             const int n
         ) {
             const int b = get_global_id(0);
+            const int base = b * n;
             const int i0 = get_global_id(1) * 4;
-            for (int i = i0; i < i0 + 4 && i < n; i++) {
+            const int i_end = min(i0 + 4, n);
+            for (int i = i0; i < i_end; i++) {
                 float2 gv = g[i];
-                int idx = b * n + i;
-                float2 a = fft_batch[idx];
-                fft_batch[idx] = (float2)(a.x * gv.x - a.y * gv.y,
-                                          a.x * gv.y + a.y * gv.x);
+                float2 a = fft_batch[base + i];
+                fft_batch[base + i] = (float2)(fma(a.x, gv.x, -(a.y * gv.y)),
+                                               fma(a.x, gv.y,   a.y * gv.x));
             }
         }
         """
@@ -1752,10 +1754,13 @@ class PoissonSolver(object):
                 'src_plane': src_plane,
                 'dst_plane': dst_plane,
                 'f_batch_dev': cla.zeros(queue, (2,) + full_shape, dtype=numpy.float32),
-                'host_pad': numpy.zeros((2,) + full_shape, dtype=numpy.float32),
+                'host_stage': numpy.empty((2,) + fshape, dtype=numpy.float32),
+                'f_stage_dev': cla.empty(queue, (2,) + fshape, dtype=numpy.float32),
                 'f_fft_batch_dev': cla.empty(queue, (2,) + fft_shape, dtype=numpy.complex64),
                 'conv_batch_dev': cla.empty(queue, (2,) + full_shape, dtype=numpy.float32),
                 'p_batch_dev': cla.empty(queue, (2,) + fshape, dtype=numpy.float32, allocator=alloc_pinned),
+                'p1_host': numpy.empty(fshape, dtype=numpy.float32),
+                'p2_host': numpy.empty(fshape, dtype=numpy.float32),
             }
             timing['cache_prepare_ms'] = (timeit.default_timer() - t0) * 1000.0
 
@@ -1763,15 +1768,26 @@ class PoissonSolver(object):
         fshape_c = c['fshape']
         full_shape_c = c['full_shape']
 
-        # H2D: write into pre-zeroed host pad buffer, single enqueue_copy (no map/unmap/scatter)
+        # H2D: compact transfer of unpadded data, then scatter to padded buffer on device
         t0 = timeit.default_timer()
-        host_pad = c['host_pad']
         H, W = fshape_c
-        f1_f32 = numpy.ascontiguousarray(f1, dtype=numpy.float32)
-        f2_f32 = numpy.ascontiguousarray(f2, dtype=numpy.float32)
-        host_pad[0, :H, :W] = f1_f32
-        host_pad[1, :H, :W] = f2_f32
-        cl.enqueue_copy(queue, c['f_batch_dev'].data, host_pad, is_blocking=False)
+        f1_f32 = f1 if (f1.flags['C_CONTIGUOUS'] and f1.dtype == numpy.float32) else numpy.ascontiguousarray(f1, dtype=numpy.float32)
+        f2_f32 = f2 if (f2.flags['C_CONTIGUOUS'] and f2.dtype == numpy.float32) else numpy.ascontiguousarray(f2, dtype=numpy.float32)
+        host_stage = c['host_stage']
+        numpy.copyto(host_stage[0], f1_f32)
+        numpy.copyto(host_stage[1], f2_f32)
+        f_stage_dev = c['f_stage_dev']
+        cl.enqueue_copy(queue, f_stage_dev.data, host_stage, is_blocking=False)
+        c['f_batch_dev'].fill(0, queue=queue)
+        scatter_batch = self._get_opencl_scatter_batch_kernel(queue)
+        scatter_batch(
+            queue, (W, H, 2), None,
+            f_stage_dev.data, c['f_batch_dev'].data,
+            numpy.int32(W), numpy.int32(H),
+            numpy.int32(full_shape_c[1]),
+            numpy.int32(H * W),
+            numpy.int32(c['dst_plane']),
+        )
         if profile:
             queue.finish()
         timing['h2d_set_ms'] = (timeit.default_timer() - t0) * 1000.0
@@ -1827,11 +1843,13 @@ class PoissonSolver(object):
             queue.finish()
         timing['crop_scale_ms'] = (timeit.default_timer() - t0) * 1000.0
 
-        # D2H: map pinned buffer (ALLOC_HOST_PTR) for read, copy out p1/p2, unmap
+        # D2H: map pinned buffer, copy into pre-allocated host arrays (zero allocation)
         t0 = timeit.default_timer()
+        p1 = c['p1_host']
+        p2 = c['p2_host']
         mapped_p = c['p_batch_dev'].map_to_host(queue=queue, flags=cl.map_flags.READ, is_blocking=True)
-        p1 = numpy.asarray(mapped_p[0], dtype=numpy.float32).copy()
-        p2 = numpy.asarray(mapped_p[1], dtype=numpy.float32).copy()
+        numpy.copyto(p1, mapped_p[0])
+        numpy.copyto(p2, mapped_p[1])
         mapped_p.base.release(queue=queue)
         timing['d2h_get_ms'] = (timeit.default_timer() - t0) * 1000.0
 
diff --git a/src/orisys/Stitcher.py b/src/orisys/Stitcher.py
index b63434f..e6e4cfd 100644
--- a/src/orisys/Stitcher.py
+++ b/src/orisys/Stitcher.py
@@ -406,7 +406,7 @@ class Stitcher:
             wx[-m:] = ramp[::-1]
         return wy[:, None] * wx[None, :]
 
-    def _blend_tile_inplace(self, blend, wsum, cov, src, weight, y0, x0):
+    def _blend_tile_inplace(self, blend, wsum, cov, src, weight, y0, x0, scratch=None):
         """Copy non-overlap tiles directly; blend only when this tile touches existing content."""
         h, w = src.shape[:2]
         H, W = blend.shape[:2]
@@ -426,6 +426,19 @@ class Stitcher:
         if not np.any(cov_patch):
             blend_patch[:] = src_patch
             wsum_patch[:] = w_patch
+        elif scratch is not None:
+            s_src, s_alpha = scratch
+            sf = s_src[:sh, :sw]
+            sa = s_alpha[:sh, :sw]
+            sf[:] = src_patch
+            np.add(wsum_patch, w_patch, out=wsum_patch)
+            np.maximum(wsum_patch, 1e-6, out=sa)
+            np.divide(w_patch, sa, out=sa)
+            np.subtract(sf, blend_patch, out=sf)
+            sf[:, :, 0] *= sa
+            sf[:, :, 1] *= sa
+            sf[:, :, 2] *= sa
+            np.add(blend_patch, sf, out=blend_patch)
         else:
             total = wsum_patch + w_patch
             alpha = w_patch / np.maximum(total, 1e-6)
@@ -436,13 +449,15 @@ class Stitcher:
         cov_patch[:] = 1
 
     def _mean_gray_bgr(self, bgr):
-        return float(cv2.cvtColor(bgr, cv2.COLOR_BGR2GRAY).mean())
+        ch = cv2.mean(bgr)
+        return (ch[0] + ch[1] + ch[2]) / 3.0
 
     def _match_bgr_brightness(self, tile_bgr, ref_mean_gray):
         delta = int(round(ref_mean_gray - self._mean_gray_bgr(tile_bgr)))
         if delta == 0:
             return tile_bgr
-        return np.clip(tile_bgr.astype(np.int16) + delta, 0, 255).astype(np.uint8)
+        lut = np.clip(np.arange(256, dtype=np.int16) + delta, 0, 255).astype(np.uint8)
+        return cv2.LUT(tile_bgr, lut)
 
     def _largest_fully_covered_rect(self, cov):
         if cov.size == 0:
@@ -519,13 +534,22 @@ class Stitcher:
         do_profile = verbose and frame_idx % print_every == 0
         t_profile0 = cv2.getTickCount() / cv2.getTickFrequency() if do_profile else 0.0
 
-        t_ref_bgr = self._tile_to_bgr(t_ref_bgr)
+        if t_ref_bgr.ndim == 2:
+            t_ref_bgr = cv2.cvtColor(t_ref_bgr, cv2.COLOR_GRAY2BGR)
         h0c, w0c = t_ref_bgr.shape[:2]
         if match_brightness:
             ref_mean_gray = self._mean_gray_bgr(t_ref_bgr)
-            tiles_bgr_c = [self._match_bgr_brightness(self._tile_to_bgr(tb), ref_mean_gray) for tb in tiles_bgr]
+            tiles_bgr_c = [
+                self._match_bgr_brightness(
+                    tb if tb.ndim == 3 else cv2.cvtColor(tb, cv2.COLOR_GRAY2BGR),
+                    ref_mean_gray,
+                )
+                for tb in tiles_bgr
+            ]
+        elif any(tb.ndim == 2 for tb in tiles_bgr):
+            tiles_bgr_c = [tb if tb.ndim == 3 else cv2.cvtColor(tb, cv2.COLOR_GRAY2BGR) for tb in tiles_bgr]
         else:
-            tiles_bgr_c = [self._tile_to_bgr(tb) for tb in tiles_bgr]
+            tiles_bgr_c = tiles_bgr
         t_after_brightness = cv2.getTickCount() / cv2.getTickFrequency() if do_profile else 0.0
 
         geom_key = self._geom_cache_key(t_ref_bgr, tiles_bgr_c, pad)
@@ -563,6 +587,20 @@ class Stitcher:
                 tx, ty = self._translation_from_H_to_ref(H_use)
                 pastes.append((int(np.floor(ty - y_min_f)), int(np.floor(tx - x_min_f))))
 
+            _dr_y0 = max(0, ref_y0)
+            _dr_x0 = max(0, ref_x0)
+            _dr_y1 = min(out_h, ref_y0 + h0c)
+            _dr_x1 = min(out_w, ref_x0 + w0c)
+            for _tb, _pa in zip(tiles_bgr_c, pastes):
+                if _pa is None:
+                    continue
+                _py, _px = _pa
+                _th, _tw = _tb.shape[:2]
+                _dr_y0 = min(_dr_y0, max(0, _py))
+                _dr_x0 = min(_dr_x0, max(0, _px))
+                _dr_y1 = max(_dr_y1, min(out_h, _py + _th))
+                _dr_x1 = max(_dr_x1, min(out_w, _px + _tw))
+
             canvas_cache["geom_key"] = geom_key
             canvas_cache["out_w"] = out_w
             canvas_cache["out_h"] = out_h
@@ -572,6 +610,7 @@ class Stitcher:
             canvas_cache["x_min_f"] = x_min_f
             canvas_cache["y_min_f"] = y_min_f
             canvas_cache["crop_rect"] = None
+            canvas_cache["dirty_rect"] = (_dr_y0, _dr_y1, _dr_x0, _dr_x1)
             canvas_cache["blend_buf"] = np.zeros((out_h, out_w, 3), dtype=np.float32)
             canvas_cache["wsum_buf"] = np.zeros((out_h, out_w), dtype=np.float32)
             canvas_cache["cov_buf"] = np.zeros((out_h, out_w), dtype=np.uint8)
@@ -588,12 +627,32 @@ class Stitcher:
         wsum = canvas_cache["wsum_buf"]
         cov = canvas_cache["cov_buf"]
         weight_maps = canvas_cache["weight_maps"]
-        blend_buf.fill(0)
-        wsum.fill(0)
-        cov.fill(0)
+        if geom_cache_hit:
+            _dr = canvas_cache.get("dirty_rect")
+            if _dr is not None:
+                _dy0, _dy1, _dx0, _dx1 = _dr
+                blend_buf[_dy0:_dy1, _dx0:_dx1] = 0
+                wsum[_dy0:_dy1, _dx0:_dx1] = 0
+                cov[_dy0:_dy1, _dx0:_dx1] = 0
+            else:
+                blend_buf.fill(0)
+                wsum.fill(0)
+                cov.fill(0)
         t_after_reset = cv2.getTickCount() / cv2.getTickFrequency() if do_profile else 0.0
 
         blend_margin = 16 if blend else 0
+        max_th = max(h0c, *(tb.shape[0] for tb in tiles_bgr_c)) if tiles_bgr_c else h0c
+        max_tw = max(w0c, *(tb.shape[1] for tb in tiles_bgr_c)) if tiles_bgr_c else w0c
+        _ss = canvas_cache.get("scratch_src")
+        _sa = canvas_cache.get("scratch_alpha")
+        if _ss is None or _ss.shape[0] < max_th or _ss.shape[1] < max_tw:
+            _ss = np.empty((max_th, max_tw, 3), dtype=np.float32)
+            canvas_cache["scratch_src"] = _ss
+        if _sa is None or _sa.shape[0] < max_th or _sa.shape[1] < max_tw:
+            _sa = np.empty((max_th, max_tw), dtype=np.float32)
+            canvas_cache["scratch_alpha"] = _sa
+        _scratch = (_ss, _sa)
+
         blended_tiles = 0
         for tile_bgr, H_use, paste in zip(tiles_bgr_c, H_use_list, pastes):
             if H_use is None or paste is None:
@@ -603,13 +662,13 @@ class Stitcher:
             wkey = (th, tw, blend_margin)
             if wkey not in weight_maps:
                 weight_maps[wkey] = self._tile_weight_map(th, tw, blend_margin)
-            self._blend_tile_inplace(blend_buf, wsum, cov, tile_bgr, weight_maps[wkey], y0, x0)
+            self._blend_tile_inplace(blend_buf, wsum, cov, tile_bgr, weight_maps[wkey], y0, x0, _scratch)
             blended_tiles += 1
 
         ref_wkey = (h0c, w0c, blend_margin)
         if ref_wkey not in weight_maps:
             weight_maps[ref_wkey] = self._tile_weight_map(h0c, w0c, blend_margin)
-        self._blend_tile_inplace(blend_buf, wsum, cov, t_ref_bgr, weight_maps[ref_wkey], ref_y0, ref_x0)
+        self._blend_tile_inplace(blend_buf, wsum, cov, t_ref_bgr, weight_maps[ref_wkey], ref_y0, ref_x0, _scratch)
         t_after_blend = cv2.getTickCount() / cv2.getTickFrequency() if do_profile else 0.0
 
         crop_rect = canvas_cache.get("crop_rect")
@@ -618,9 +677,16 @@ class Stitcher:
             canvas_cache["crop_rect"] = crop_rect
         if crop_rect is not None:
             y0, y1, x0, x1 = crop_rect
-            stitched = np.clip(blend_buf[y0:y1, x0:x1], 0, 255).astype(np.uint8)
+            roi = blend_buf[y0:y1, x0:x1]
         else:
-            stitched = np.clip(blend_buf, 0, 255).astype(np.uint8)
+            roi = blend_buf
+        rh, rw = roi.shape[:2]
+        stitch_buf = canvas_cache.get("stitch_buf")
+        if stitch_buf is None or stitch_buf.shape[0] != rh or stitch_buf.shape[1] != rw:
+            stitch_buf = np.empty((rh, rw, 3), dtype=np.uint8)
+            canvas_cache["stitch_buf"] = stitch_buf
+        cv2.convertScaleAbs(roi, dst=stitch_buf)
+        stitched = stitch_buf
         t_after_finalize = cv2.getTickCount() / cv2.getTickFrequency() if do_profile else 0.0
 
         if do_profile:
diff --git a/src/orisys/sensor.py b/src/orisys/sensor.py
index 8f415f6..ba54a9a 100644
--- a/src/orisys/sensor.py
+++ b/src/orisys/sensor.py
@@ -60,7 +60,8 @@ class Sensor:
     # =========================================================================
     
     def __init__(self, vid_src, cuda=True, isstitch=True, verbose=True, debug=False,
-                 config_name="prototype", cal_path=None, override=False, box=None, mirror=True, backend='auto'):
+                 config_name="prototype", cal_path=None, override=False, box=None, mirror=True,
+                 backend='auto'):
         """
         Initialize Sensor with video source and configuration.
         
@@ -78,6 +79,8 @@ class Sensor:
                 "cupy": Use CuPy for GPU acceleration
                 "opencl": Use OpenCL for GPU acceleration
                 "cpu": Use CPU for computation
+            hhd_downsample: Downsample factor for flow field before HHD (default: 1, no downsampling).
+                Factor of 2 halves each dimension, giving ~4x speedup in Poisson solve.
         """
 
         self.CUPY_AVAILABLE = check_cupy_available()
@@ -147,7 +150,12 @@ class Sensor:
         self.__tracker = Tracker(adaptive=True, cuda=cuda, verbose=self.verbose, config_tracker=self.__config_tracker)
         self.__contact3d = Contact3D(config_contact=self.__config_contact, verbose=self.verbose)
         self.__hhd = None
-        self.poisson_backend = backend
+        # cfg_ds = self.__config_tracker.get("hhd_downsample", None)
+        # self._hhd_downsample = max(1, int(cfg_ds if cfg_ds is not None else hhd_downsample))
+        self._hhd_downsample = self.__config_tracker.get("hhd_downsample", None)
+        self._hhd_grid_shape = None
+        cfg_backend = self.__config_tracker.get("backend", None)
+        self.poisson_backend = cfg_backend if cfg_backend is not None else backend
         self._use_pyav = False
         self._av_container = None
         self._av_stream = None
@@ -287,28 +295,36 @@ class Sensor:
         flow = self.flow
         if flow is None:
             raise Exception("Warning: Compute flow first")
-        
-        # Initialize HHD solver
-        if self.__hhd is None:
-            self.__hhd = nHHD(cuda=self.cuda, grid=flow.shape[:2], spacings=(18, 18), poisson_backend=self.poisson_backend, debug=self.debug)
 
-        # Force computation
-        # When OpenCL backend is selected, use host force path to ensure
-        # Poisson solve goes through PoissonSolver.solve_opencl().
-        _, forces = self.__hhd.get_force(flow)
-        # if self.poisson_backend == 'opencl':
-        #     _, forces = self.__hhd.get_force_host(flow)
-        # elif self.cuda and self.CUPY_AVAILABLE:
-        #     _, forces = self.__hhd.get_force_gpu_v2(flow)
-        # else:
-        #     _, forces = self.__hhd.get_force_host(flow)
+        orig_h, orig_w = flow.shape[:2]
+        k = self._hhd_downsample
+
+        if k > 1:
+            ds_w = max(1, orig_w // k)
+            ds_h = max(1, orig_h // k)
+            flow_ds = cv2.resize(flow, (ds_w, ds_h), interpolation=cv2.INTER_AREA)
+        else:
+            flow_ds = flow
+
+        ds_shape = flow_ds.shape[:2]
+        if self.__hhd is None or self._hhd_grid_shape != ds_shape:
+            spacings = (18 * k, 18 * k)
+            self.__hhd = nHHD(cuda=self.cuda, grid=ds_shape, spacings=spacings, poisson_backend=self.poisson_backend, debug=self.debug)
+            self._hhd_grid_shape = ds_shape
+
+        _, forces = self.__hhd.get_force(flow_ds)
         self.fnormal, self.fshearx, self.fsheary, self.ftor = forces
-        
-        # Get components
-        self.vnormal = self.__hhd.get_normal()
-        self.vshear = self.__hhd.get_shear()
 
-        # Update state
+        vnormal_ds = self.__hhd.get_normal()
+        vshear_ds = self.__hhd.get_shear()
+
+        if k > 1:
+            self.vnormal = cv2.resize(vnormal_ds, (orig_w, orig_h), interpolation=cv2.INTER_LINEAR)
+            self.vshear = cv2.resize(vshear_ds, (orig_w, orig_h), interpolation=cv2.INTER_LINEAR)
+        else:
+            self.vnormal = vnormal_ds
+            self.vshear = vshear_ds
+
         self.compute_deformation_valid = True
 
         return True
diff --git a/src/orisys/tracker.py b/src/orisys/tracker.py
index b6db0f0..5a8a67a 100644
--- a/src/orisys/tracker.py
+++ b/src/orisys/tracker.py
@@ -63,8 +63,8 @@ class Tracker:
             # self.__tracking_inst.setUseSpatialPropagation(False)
             # self.__tracking_inst.setUseMeanNormalization(False)
             # self.__tracking_inst.setGradientDescentIterations(16) #16
-            self.__tracking_inst.setPatchSize(8) # 8
-            self.__tracking_inst.setPatchStride(1) # 3
+            # self.__tracking_inst.setPatchSize(12) # 8
+            # self.__tracking_inst.setPatchStride(5) # 3
             # self.__tracking_inst.setVariationalRefinementIterations(5)#5
 
             # else:
diff --git a/test_downsample.py b/test_downsample.py
new file mode 100644
index 0000000..d9e1f77
--- /dev/null
+++ b/test_downsample.py
@@ -0,0 +1,252 @@
+"""
+HHD Downsample Factor Comparison
+
+Reads a video file and runs the full pipeline at multiple hhd_downsample factors.
+Computes per-frame HHD timing and vnormal quality metrics against the baseline (factor=1),
+then plots the results.
+
+Usage:
+    python test_downsample.py -v jx01_hex1.5_2.mp4 -c ./config/prototype_of3_720p.json -f 1,2,3,4 -n 50 -w 30
+"""
+import os
+import sys
+
+os.environ["OPENCV_VIDEOIO_MSMF_ENABLE_HW_TRANSFORMS"] = "0"
+
+_local_src = os.path.join(os.path.dirname(__file__), "src")
+if _local_src not in sys.path:
+    sys.path.insert(0, _local_src)
+
+import cv2
+import orisys
+import numpy as np
+import time
+import argparse
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+
+
+def run_sensor(vid_src, config, cal, downsample, n_warmup, n_measure):
+    """Run sensor for n_warmup + n_measure frames. Capture timing and vnormal for every measurement frame."""
+    try:
+        vid = int(vid_src)
+    except ValueError:
+        vid = vid_src
+
+    sensor = orisys.Sensor(
+        vid, config_name=config, cal_path=cal,
+        isstitch=True, cuda=False, verbose=False,
+        hhd_downsample=downsample,
+    )
+
+    hhd_times = []
+    vnormals = []
+    frame_idx = 0
+    n_total = n_warmup + n_measure
+
+    while frame_idx < n_total:
+        try:
+            sensor.get_img()
+        except Exception:
+            break
+
+        sensor.compute_flow()
+
+        t0 = time.perf_counter()
+        sensor.compute_hhd()
+        t_hhd = time.perf_counter() - t0
+
+        if frame_idx >= n_warmup:
+            hhd_times.append(t_hhd)
+            vn = sensor.vnormal.copy() if sensor.vnormal is not None else None
+            vnormals.append(vn)
+
+        frame_idx += 1
+
+    sensor.disconnect()
+    return {"hhd_times": hhd_times, "vnormals": vnormals}
+
+
+def per_frame_rmse(a, b):
+    if a is None or b is None:
+        return float("nan")
+    diff = a.astype(np.float64) - b.astype(np.float64)
+    return float(np.sqrt(np.mean(diff ** 2)))
+
+
+def per_frame_max_abs_err(a, b):
+    if a is None or b is None:
+        return float("nan")
+    return float(np.max(np.abs(a.astype(np.float64) - b.astype(np.float64))))
+
+
+def per_frame_cosine_sim(a, b):
+    if a is None or b is None:
+        return float("nan")
+    a64 = a.astype(np.float64)
+    b64 = b.astype(np.float64)
+    dot = np.sum(a64 * b64, axis=-1)
+    mag_a = np.linalg.norm(a64, axis=-1)
+    mag_b = np.linalg.norm(b64, axis=-1)
+    denom = mag_a * mag_b
+    valid = denom > 1e-12
+    if not np.any(valid):
+        return float("nan")
+    return float(np.mean(dot[valid] / denom[valid]))
+
+
+def main():
+    parser = argparse.ArgumentParser(description="Compare HHD downsample factors")
+    parser.add_argument("--video", "-v", type=str, default="2", help="Camera ID or video file path")
+    parser.add_argument("--config", "-c", type=str, default="./config/prototype_of3_720p.json")
+    parser.add_argument("--cal", "-cal", type=str, default="./config/ddjx01.npy")
+    parser.add_argument("--factors", "-f", type=str, default="1,2,3", help="Comma-separated downsample factors")
+    parser.add_argument("--warmup", "-w", type=int, default=5, help="Warmup frames (skipped)")
+    parser.add_argument("--frames", "-n", type=int, default=30, help="Measurement frames")
+    parser.add_argument("--output", "-o", type=str, default="downsample_comparison.png", help="Output plot filename")
+    args = parser.parse_args()
+
+    factors = [int(f) for f in args.factors.split(",")]
+
+    print(f"Video: {args.video}")
+    print(f"Config: {args.config}")
+    print(f"Factors: {factors}")
+    print(f"Warmup: {args.warmup} frames, Measure: {args.frames} frames")
+    print()
+
+    # --- collect data ---
+    results = {}
+    for factor in factors:
+        print(f"--- Running factor {factor} ---")
+        results[factor] = run_sensor(
+            args.video, args.config, args.cal,
+            factor, args.warmup, args.frames,
+        )
+        n = len(results[factor]["hhd_times"])
+        if n > 0:
+            arr = np.array(results[factor]["hhd_times"]) * 1000
+            print(f"  {n} frames: hhd mean={arr.mean():.2f}ms, median={np.median(arr):.2f}ms, "
+                  f"min={arr.min():.2f}ms, max={arr.max():.2f}ms")
+        else:
+            print("  No frames captured")
+        print()
+
+    # --- timing summary table ---
+    print("=" * 72)
+    print(f"{'Factor':<10} {'Mean HHD (ms)':<16} {'Median (ms)':<14} {'Speedup':<10}")
+    print("-" * 72)
+    baseline_median = None
+    for factor in factors:
+        arr = np.array(results[factor]["hhd_times"]) * 1000
+        if len(arr) == 0:
+            continue
+        med = float(np.median(arr))
+        if baseline_median is None:
+            baseline_median = med
+        speedup = baseline_median / med if med > 0 else 0
+        print(f"{factor:<10} {arr.mean():<16.2f} {med:<14.2f} {speedup:<10.2f}x")
+    print("=" * 72)
+    print()
+
+    # --- compute per-frame quality metrics vs baseline (factor=1) ---
+    has_baseline = 1 in factors and len(factors) >= 2
+    metrics = {}  # factor -> {"rmse": [...], "max_err": [...], "cosine": [...]}
+    if has_baseline:
+        base_vn = results[1]["vnormals"]
+        n_base = len(base_vn)
+        for factor in factors:
+            if factor == 1:
+                continue
+            test_vn = results[factor]["vnormals"]
+            n_common = min(n_base, len(test_vn))
+            rmse_arr, maxerr_arr, cos_arr = [], [], []
+            for i in range(n_common):
+                rmse_arr.append(per_frame_rmse(base_vn[i], test_vn[i]))
+                maxerr_arr.append(per_frame_max_abs_err(base_vn[i], test_vn[i]))
+                cos_arr.append(per_frame_cosine_sim(base_vn[i], test_vn[i]))
+            metrics[factor] = {
+                "rmse": np.array(rmse_arr),
+                "max_err": np.array(maxerr_arr),
+                "cosine": np.array(cos_arr),
+            }
+
+        # print per-factor summary
+        print(f"{'Factor':<10} {'Mean RMSE':<14} {'Mean MaxErr':<14} {'Mean CosSim':<14}")
+        print("-" * 52)
+        for factor in sorted(metrics):
+            m = metrics[factor]
+            print(f"{factor:<10} {np.nanmean(m['rmse']):<14.6f} {np.nanmean(m['max_err']):<14.4f} "
+                  f"{np.nanmean(m['cosine']):<14.6f}")
+        print("=" * 52)
+        print()
+
+    # --- plot ---
+    n_quality_factors = len(metrics)
+    n_plots = 1 + (3 if has_baseline else 0)  # timing + (rmse, max_err, cosine)
+    fig, axes = plt.subplots(n_plots, 1, figsize=(12, 3.5 * n_plots), sharex=True)
+    if n_plots == 1:
+        axes = [axes]
+
+    colors = plt.cm.tab10.colors
+    frame_nums = None
+
+    # 1) HHD timing per frame
+    ax = axes[0]
+    for i, factor in enumerate(factors):
+        arr = np.array(results[factor]["hhd_times"]) * 1000
+        frames = np.arange(len(arr))
+        if frame_nums is None:
+            frame_nums = frames
+        ax.plot(frames, arr, label=f"ds={factor}", color=colors[i % len(colors)], linewidth=1)
+    ax.set_ylabel("HHD time (ms)")
+    ax.set_title("Per-frame HHD computation time")
+    ax.legend(loc="upper right")
+    ax.grid(True, alpha=0.3)
+
+    if has_baseline:
+        # 2) RMSE
+        ax = axes[1]
+        for i, factor in enumerate(sorted(metrics)):
+            m = metrics[factor]
+            frames = np.arange(len(m["rmse"]))
+            ci = factors.index(factor) if factor in factors else i
+            ax.plot(frames, m["rmse"], label=f"ds={factor}", color=colors[ci % len(colors)], linewidth=1)
+        ax.set_ylabel("RMSE (vnormal)")
+        ax.set_title("Per-frame RMSE vs baseline (ds=1)")
+        ax.legend(loc="upper right")
+        ax.grid(True, alpha=0.3)
+
+        # 3) Max absolute error
+        ax = axes[2]
+        for i, factor in enumerate(sorted(metrics)):
+            m = metrics[factor]
+            frames = np.arange(len(m["max_err"]))
+            ci = factors.index(factor) if factor in factors else i
+            ax.plot(frames, m["max_err"], label=f"ds={factor}", color=colors[ci % len(colors)], linewidth=1)
+        ax.set_ylabel("Max Abs Error (vnormal)")
+        ax.set_title("Per-frame max absolute error vs baseline (ds=1)")
+        ax.legend(loc="upper right")
+        ax.grid(True, alpha=0.3)
+
+        # 4) Cosine similarity
+        ax = axes[3]
+        for i, factor in enumerate(sorted(metrics)):
+            m = metrics[factor]
+            frames = np.arange(len(m["cosine"]))
+            ci = factors.index(factor) if factor in factors else i
+            ax.plot(frames, m["cosine"], label=f"ds={factor}", color=colors[ci % len(colors)], linewidth=1)
+        ax.set_ylabel("Cosine Similarity")
+        ax.set_title("Per-frame cosine similarity vs baseline (ds=1)")
+        ax.legend(loc="upper right")
+        ax.grid(True, alpha=0.3)
+        ax.set_ylim(bottom=min(0.999, ax.get_ylim()[0]))
+
+    axes[-1].set_xlabel("Frame")
+    fig.tight_layout()
+    fig.savefig(args.output, dpi=150)
+    print(f"Plot saved to: {args.output}")
+
+
+if __name__ == "__main__":
+    main()