Image wave refactor

.vscode/c_cpp_properties.json

@@ -3,12 +3,14 @@
         {
             "name": "Linux",
             "includePath": [
-                "${workspaceFolder}/**"
+                "${workspaceFolder}/**",
+                "${workspaceFolder}/src/vstarstack/library/fine_movement/libimagedeform/include",
+                "/usr/include/python3.11"
             ],
             "defines": [],
             "compilerPath": "/usr/bin/gcc",
-            "cStandard": "gnu17",
-            "cppStandard": "gnu++17",
+            "cStandard": "c17",
+            "cppStandard": "c++17",
             "intelliSenseMode": "linux-gcc-x64"
         }
     ],

.vscode/settings.json

@@ -4,7 +4,13 @@
         "random": "c",
         "math.h": "c",
         "unistd.h": "c",
-        "sphere.h": "c"
+        "sphere.h": "c",
+        "image_grid.h": "c",
+        "image_deform.h": "c",
+        "image_deform_gc.h": "c",
+        "stdlib.h": "c",
+        "structmember.h": "c",
+        "imagedeform_gc.h": "c"
     },
     "editor.tabSize": 4,
     "editor.insertSpaces": true,

setup.py

@@ -35,14 +35,26 @@
                             np.get_include(),
                         ])
 
-image_wave = Extension(name="vstarstack.library.fine_shift.image_wave",
-       sources=["src/vstarstack/library/fine_shift/image_wave.c",
-                "src/vstarstack/library/fine_shift/image_wave_interpolation.c",
-                "src/vstarstack/library/fine_shift/image_wave_targets.c",
-                "src/vstarstack/library/fine_shift/image_wave_correlation.c",
-                "src/vstarstack/library/fine_shift/image_wave_module.c",
-                "src/vstarstack/library/fine_shift/image_wave_image.c",
-                ], include_dirs=[np.get_include()])
+libimagedeform_root = "src/vstarstack/library/fine_movement/libimagedeform"
+libimagedeform_headers = [libimagedeform_root + "/include"]
+libimagedeform_sources = [libimagedeform_root + "/src/interpolation.c",
+                          libimagedeform_root + "/src/image_grid.c",
+                          libimagedeform_root + "/src/image_deform.c",
+                          libimagedeform_root + "/src/image_deform_gc.c",
+                          libimagedeform_root + "/src/image_deform_lc.c",
+                        ]
+
+imagedeform_root = "src/vstarstack/library/fine_movement/module"
+imagedeform_sources = [imagedeform_root + "/imagegrid.c",
+                       imagedeform_root + "/imagedeform.c",
+                       imagedeform_root + "/imagedeform_gc.c",
+                       imagedeform_root + "/imagedeform_lc.c",
+                       imagedeform_root + "/imagedeform_module.c",
+                        ]
+
+image_deform = Extension(name="vstarstack.library.fine_movement.module",
+       sources=imagedeform_sources+libimagedeform_sources,
+       include_dirs=[np.get_include()]+libimagedeform_headers)
 
 root = os.path.join(os.path.abspath(os.path.dirname(__file__)), "src")
 result = [os.path.join(dp, f) for dp, dn, filenames in os.walk(root)
@@ -63,7 +75,8 @@
        packages=packages,
        ext_modules = [projection,
                       movements,
-                      image_wave],
+                      image_deform,
+                    ],
        install_requires = [
               'numpy',
               'astropy',

src/vstarstack/library/fine_movement/aligner.py

@@ -0,0 +1,184 @@
+#
+# Copyright (c) 2022-2024 Vladislav Tsendrovskii
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, version 3 of the License.
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
+# See the GNU General Public License for more details.
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see <https://www.gnu.org/licenses/>.
+#
+
+import numpy as np
+import scipy
+
+from vstarstack.library.fine_movement.module import ImageGrid
+from vstarstack.library.fine_movement.module import ImageDeform
+from vstarstack.library.fine_movement.module import ImageDeformLC
+from vstarstack.library.fine_movement.module import ImageDeformGC
+
+from vstarstack.library.data import DataFrame
+
+def _cluster_average(cluster):
+    xs = [cluster[name]["x"] for name in cluster]
+    ys = [cluster[name]["y"] for name in cluster]
+    av = {
+        "x" : sum(xs) / len(xs),
+        "y" : sum(ys) / len(ys),
+    }
+    return av
+
+class Aligner:
+    """Alignment applier"""
+
+    def __init__(self, image_w : int, image_h : int, shift_array : np.ndarray):
+        self.image_w = image_w
+        self.image_h = image_h
+        self.grid_h = shift_array.shape[0]
+        self.grid_w = shift_array.shape[1]
+        self.deform = ImageDeform(self.image_w, self.image_h, self.grid_w, self.grid_h)
+        self.deform.fill(shift_array=shift_array)
+
+    def apply_alignment(self, dataframe : DataFrame, subpixels : int):
+        """Apply alignment descriptor to file"""
+        for channel in dataframe.get_channels():
+            image, opts = dataframe.get_channel(channel)
+            if opts["encoded"]:
+                continue
+            w = image.shape[1]
+            h = image.shape[0]
+            grid = ImageGrid(w, h)
+            grid.fill(image.astype('double'))
+            fixed_grid = self.deform.apply_image(image=grid, subpixels=subpixels)
+            fixed = fixed_grid.content()
+            fixed[np.where(np.isnan(fixed))] = 0
+            dataframe.replace_channel(fixed, channel)
+        return dataframe
+
+    def serialize(self):
+        """Serialize image deform"""
+        nitems = self.grid_h * self.grid_w * 2
+        shift_array = np.reshape(self.deform.content(), (nitems,))
+        return {
+            "grid_w" : self.grid_w,
+            "grid_h" : self.grid_h,
+            "image_w" : self.image_w,
+            "image_h" : self.image_h,
+            "array" : list(shift_array),
+        }
+
+    @staticmethod
+    def deserialize(description : dict):
+        """Deserialize image deform"""
+        grid_h = description["grid_h"]
+        grid_w = description["grid_w"]
+        shift_array = np.reshape(description["array"], (grid_h, grid_w, 2))
+        return Aligner(description["image_w"], description["image_h"], shift_array)
+
+class ClusterAlignerBuilder:
+    def __init__(self, image_w, image_h, grid_w, grid_h, spk, num_steps, min_points, dh):
+        self.image_w = image_w
+        self.image_h = image_h
+        self.grid_w = grid_w
+        self.grid_h = grid_h
+        self.spk = spk
+        self.num_steps = num_steps
+        self.min_points = min_points
+        self.dh = dh
+        self.correlator = ImageDeformGC(self.image_w, self.image_h,
+                                        self.grid_w, self.grid_h,
+                                        self.spk)
+
+    def find_alignment(self, name : str, clusters : list) -> Aligner:
+        """Find alignment of image `name` using clusters"""
+        expected_points = []
+        actual_points = []
+        for cluster in clusters:
+            if name not in cluster:
+                continue
+            average = _cluster_average(cluster)
+
+            # we need reverse transformation
+            x = average["x"]
+            y = average["y"]
+            expected_points.append((x, y))
+            x = cluster[name]["x"]
+            y = cluster[name]["y"]
+            actual_points.append((x, y))
+
+        print(f"\tusing {len(expected_points)} points")
+        if len(expected_points) < self.min_points:
+            print("\tskip - too low points")
+            return None
+
+        deform = self.correlator.find(points=actual_points,
+                                      expected_points=expected_points,
+                                      dh=self.dh,
+                                      Nsteps=self.num_steps)
+        shift_array = deform.content()
+        return Aligner(self.image_w, self.image_h, shift_array)
+
+    def find_all_alignments(self, clusters : list) -> dict:
+        """Build alignment descriptor using clusters"""
+        names = []
+        for cluster in clusters:
+            names += cluster.keys()
+        names = set(names)
+        deforms = {}
+        for name in names:
+            deform = self.find_alignment(name, clusters)
+            if deform is not None:
+                deforms[name] = deform
+        return deforms
+
+class CorrelationAlignedBuilder:
+
+    def __init__(self,
+                 image_w : int,
+                 image_h : int,
+                 pixels : int,
+                 radius : int,
+                 maximal_shift : float,
+                 subpixels : int):
+        self.image_w = image_w
+        self.image_h = image_h
+        self.pixels = pixels
+        self.correlator = ImageDeformLC(self.image_w, self.image_h, self.pixels)
+
+        self.radius = radius
+        self.max_shift = maximal_shift
+        self.subpixels = subpixels
+
+    def find_alignment(self, image : np.ndarray,
+                       image_ref : np.ndarray,
+                       pre_align : Aligner | None,
+                       pre_align_ref : Aligner | None,
+                       smooth : int | None) -> Aligner:
+        """Build alignment descriptor of image using correlations"""
+        if image.shape != image_ref.shape:
+            return None
+        w = image.shape[1]
+        h = image.shape[0]
+
+        grid = ImageGrid(w, h)
+        grid.fill(image)
+        grid_ref = ImageGrid(w, h)
+        grid_ref.fill(image_ref)
+
+        deform_img = pre_align.deform if pre_align is not None else None
+        deform_ref = pre_align_ref.deform if pre_align_ref is not None else None
+        deform = self.correlator.find(grid, deform_img,
+                                      grid_ref, deform_ref,
+                                      self.radius,
+                                      self.max_shift,
+                                      self.subpixels)
+        print(f"smooth = {smooth}")
+        if smooth is not None:
+            data = deform.content()
+            data = scipy.ndimage.gaussian_filter(data, sigma=smooth, axes=(0,1))
+            deform.fill(data)
+        shift_array = deform.content()
+        return Aligner(self.image_w, self.image_h, shift_array)

src/vstarstack/library/fine_movement/libimagedeform/CMakeLists.txt

@@ -0,0 +1,59 @@
+cmake_minimum_required(VERSION 3.14)
+
+project(imagedeform)
+
+# Library sources
+add_library(imagedeform STATIC src/interpolation.c
+                               src/image_grid.c
+                               src/image_deform.c
+                               src/image_deform_gc.c
+                               src/image_deform_lc.c)
+
+set(HEADERS "include/interpolation.h;"
+            "include/image_grid.h;"
+            "include/image_deform.h;"
+            "include/image_deform_gc.h;"
+            "include/image_deform_lc.h")
+
+target_include_directories(imagedeform PUBLIC include)
+set_target_properties(imagedeform PROPERTIES PUBLIC_HEADER "${HEADERS}")
+
+# Fail on warnings
+if(MSVC)
+  target_compile_options(imagedeform PRIVATE /W4 /WX)
+else()
+  target_compile_options(imagedeform PRIVATE -Wall -Wextra -Wpedantic -Werror)
+endif()
+
+install(TARGETS imagedeform
+        PUBLIC_HEADER DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/libimagedeform)
+
+# Tests
+if (BUILD_TESTS)
+    set(INSTALL_GTEST OFF)
+    include(FetchContent)
+    FetchContent_Declare(
+        googletest
+        URL https://github.com/google/googletest/archive/03597a01ee50ed33e9dfd640b249b4be3799d395.zip
+    )
+
+    # For Windows: Prevent overriding the parent project's compiler/linker settings
+    set(gtest_force_shared_crt ON CACHE BOOL "" FORCE)
+    FetchContent_MakeAvailable(googletest)
+
+    enable_testing()
+    add_executable(
+        deform_test
+        tests/deform_test.cc
+        tests/deform_test_gc.cc
+    )
+
+    target_link_libraries(
+        deform_test
+        GTest::gtest_main
+        imagedeform
+    )
+
+    include(GoogleTest)
+    gtest_discover_tests(deform_test)
+endif()