From 84c5572c5f81bba8cfd73e8b64cadbd401d489c8 Mon Sep 17 00:00:00 2001
From: Tomas Farago <sensej007@email.cz>
Date: Thu, 9 Sep 2021 14:41:01 +0200
Subject: stripe-filter: add vertical sigma parameter

Which allows filtering of rings which are not perfectly straight or
which are broader. vertical-sigma parameter roughly specifies the
inclination from "perfectly vertical" stripes.
---
 src/kernels/filter.cl         | 17 ++++++++++++-----
 src/ufo-filter-stripes-task.c | 41 +++++++++++++++++++++++++++++------------
 2 files changed, 41 insertions(+), 17 deletions(-)

(limited to 'src')

diff --git a/src/kernels/filter.cl b/src/kernels/filter.cl
index ddc91f6..186b2e6 100644
--- a/src/kernels/filter.cl
+++ b/src/kernels/filter.cl
@@ -31,7 +31,8 @@ filter (global float *input,
 kernel void
 stripe_filter (global float *input,
                global float *output,
-               const float sigma)
+               const float horizontal_sigma,
+               const float vertical_sigma)
 {
     const int idx = get_global_id(0);
     const int idy = get_global_id(1);
@@ -41,11 +42,17 @@ stripe_filter (global float *input,
     /* Swap frequencies for interleaved complex input */
     /* No need to negate the second half of frequencies for Gaussian */
     const float x = idx >= width >> 1 ? (width >> 1) - (idx >> 1) : idx >> 1;
-    const float weight = exp (- x * x / (2 * sigma * sigma));
+    const float x_weight = exp (- x * x / (2 * horizontal_sigma * horizontal_sigma));
 
-    if (idy == 0) {
-        output[index] = input[index] * weight;
+    if (vertical_sigma == 0.0f) {
+        if (idy == 0) {
+            output[index] = input[index] * x_weight;
+        } else {
+            output[index] = input[index];
+        }
     } else {
-        output[index] = input[index];
+        const float y = idy >= height >> 1 ? height - idy : idy;
+        const float y_weight = exp (- y * y / (2 * vertical_sigma * vertical_sigma));
+        output[index] = input[index] * (1 - y_weight * (1 - x_weight));
     }
 }
diff --git a/src/ufo-filter-stripes-task.c b/src/ufo-filter-stripes-task.c
index bd56b7e..bb87498 100644
--- a/src/ufo-filter-stripes-task.c
+++ b/src/ufo-filter-stripes-task.c
@@ -29,7 +29,8 @@
 
 
 struct _UfoFilterStripesTaskPrivate {
-    gfloat sigma;
+    gfloat horizontal_sigma;
+    gfloat vertical_sigma;
     cl_kernel kernel;
 };
 
@@ -43,7 +44,8 @@ G_DEFINE_TYPE_WITH_CODE (UfoFilterStripesTask, ufo_filter_stripes_task, UFO_TYPE
 
 enum {
     PROP_0,
-    PROP_SIGMA,
+    PROP_HORIZONTAL_SIGMA,
+    PROP_VERTICAL_SIGMA,
     N_PROPERTIES
 };
 
@@ -76,7 +78,8 @@ ufo_filter_stripes_task_process (UfoTask *task,
 
     UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (priv->kernel, 0, sizeof (cl_mem), &in_mem));
     UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (priv->kernel, 1, sizeof (cl_mem), &out_mem));
-    UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (priv->kernel, 2, sizeof (cl_float), &priv->sigma));
+    UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (priv->kernel, 2, sizeof (cl_float), &priv->horizontal_sigma));
+    UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (priv->kernel, 3, sizeof (cl_float), &priv->vertical_sigma));
 
     profiler = ufo_task_node_get_profiler (UFO_TASK_NODE (task));
     ufo_profiler_call (profiler, cmd_queue, priv->kernel, 2, requisition->dims, NULL);
@@ -163,8 +166,11 @@ ufo_filter_stripes_task_set_property (GObject *object,
     UfoFilterStripesTaskPrivate *priv = UFO_FILTER_STRIPES_TASK_GET_PRIVATE (object);
 
     switch (property_id) {
-        case PROP_SIGMA:
-            priv->sigma = g_value_get_float (value);
+        case PROP_HORIZONTAL_SIGMA:
+            priv->horizontal_sigma = g_value_get_float (value);
+            break;
+        case PROP_VERTICAL_SIGMA:
+            priv->vertical_sigma = g_value_get_float (value);
             break;
         default:
             G_OBJECT_WARN_INVALID_PROPERTY_ID(object, property_id, pspec);
@@ -181,8 +187,11 @@ ufo_filter_stripes_task_get_property (GObject *object,
     UfoFilterStripesTaskPrivate *priv = UFO_FILTER_STRIPES_TASK_GET_PRIVATE (object);
 
     switch (property_id) {
-        case PROP_SIGMA:
-            g_value_set_float (value, priv->sigma);
+        case PROP_HORIZONTAL_SIGMA:
+            g_value_set_float (value, priv->horizontal_sigma);
+            break;
+        case PROP_VERTICAL_SIGMA:
+            g_value_set_float (value, priv->vertical_sigma);
             break;
         default:
             G_OBJECT_WARN_INVALID_PROPERTY_ID(object, property_id, pspec);
@@ -201,13 +210,20 @@ ufo_filter_stripes_task_class_init (UfoFilterStripesTaskClass *klass)
     oclass->set_property = ufo_filter_stripes_task_set_property;
     oclass->get_property = ufo_filter_stripes_task_get_property;
 
-    properties[PROP_SIGMA] =
-        g_param_spec_float ("sigma",
-            "Sigma of the gaussian window",
-            "Sigma of the gaussian window",
+    properties[PROP_HORIZONTAL_SIGMA] =
+        g_param_spec_float ("horizontal-sigma",
+            "Sigma of the gaussian window in the horizontal direction",
+            "Sigma of the gaussian window in the horizontal direction",
             0.0, G_MAXFLOAT, 1e-7,
             G_PARAM_READWRITE);
 
+    properties[PROP_VERTICAL_SIGMA] =
+        g_param_spec_float ("vertical-sigma",
+            "Sigma of the gaussian window in the vertical direction",
+            "Sigma of the gaussian window in the vertical direction",
+            0.0, G_MAXFLOAT, 0.0f,
+            G_PARAM_READWRITE);
+
     for (guint i = PROP_0 + 1; i < N_PROPERTIES; i++)
         g_object_class_install_property (oclass, i, properties[i]);
 
@@ -220,5 +236,6 @@ ufo_filter_stripes_task_init (UfoFilterStripesTask *self)
     UfoFilterStripesTaskPrivate *priv;
     self->priv = priv = UFO_FILTER_STRIPES_TASK_GET_PRIVATE (self);
     priv->kernel = NULL;
-    priv->sigma = 1e-7;
+    priv->horizontal_sigma = 1e-7;
+    priv->vertical_sigma = 0.0f;
 }
-- 
cgit v1.2.1