Add +proj=set operation to set component(s) of a coordinate to a fixed value

Fixes #1846

4 files changed, 78 insertions, 0 deletions

diff --git a/src/Makefile.am b/src/Makefile.am
index 8ab30a33..de35a754 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -174,6 +174,7 @@ libproj_la_SOURCES = \
 	conversions/geoc.cpp \
 	conversions/geocent.cpp \
 	conversions/noop.cpp \
+	conversions/set.cpp \
 	conversions/unitconvert.cpp \
 	\
 	transformations/affine.cpp \
diff --git a/src/conversions/set.cpp b/src/conversions/set.cpp
new file mode 100644
index 00000000..7628bf4f
--- /dev/null
+++ b/src/conversions/set.cpp
@@ -0,0 +1,75 @@
+#define PJ_LIB__
+
+#include "proj_internal.h"
+#include <errno.h>
+
+PROJ_HEAD(set, "Set coordinate value");
+
+/* Projection specific elements for the PJ object */
+namespace { // anonymous namespace
+struct Set {
+    bool v1;
+    bool v2;
+    bool v3;
+    bool v4;
+    double v1_val;
+    double v2_val;
+    double v3_val;
+    double v4_val;
+};
+} // anonymous namespace
+
+static PJ_COORD set_fwd_inv(PJ_COORD point, PJ *P) {
+
+    struct Set *set = static_cast<struct Set*>(P->opaque);
+
+    if (set->v1)
+        point.v[0] = set->v1_val;
+    if (set->v2)
+        point.v[1] = set->v2_val;
+    if (set->v3)
+        point.v[2] = set->v3_val;
+    if (set->v4)
+        point.v[3] = set->v4_val;
+
+    return point;
+}
+
+PJ *OPERATION(set, 0) {
+    P->inv4d = set_fwd_inv;
+    P->fwd4d = set_fwd_inv;
+
+    auto set = static_cast<struct Set*>(pj_calloc (1, sizeof(struct Set)));
+    P->opaque = set;
+    if (nullptr==P->opaque)
+        return pj_default_destructor(P, ENOMEM);
+
+    if (pj_param_exists(P->params, "v_1"))
+    {
+        set->v1 = true;
+        set->v1_val = pj_param (P->ctx, P->params, "dv_1").f;
+    }
+
+    if (pj_param_exists(P->params, "v_2"))
+    {
+        set->v2 = true;
+        set->v2_val = pj_param (P->ctx, P->params, "dv_2").f;
+    }
+
+    if (pj_param_exists(P->params, "v_3"))
+    {
+        set->v3 = true;
+        set->v3_val = pj_param (P->ctx, P->params, "dv_3").f;
+    }
+
+    if (pj_param_exists(P->params, "v_4"))
+    {
+        set->v4 = true;
+        set->v4_val = pj_param (P->ctx, P->params, "dv_4").f;
+    }
+
+    P->left  = PJ_IO_UNITS_WHATEVER;
+    P->right = PJ_IO_UNITS_WHATEVER;
+
+    return P;
+}
diff --git a/src/lib_proj.cmake b/src/lib_proj.cmake
index 704ece3d..d89cfade 100644
--- a/src/lib_proj.cmake
+++ b/src/lib_proj.cmake
@@ -187,6 +187,7 @@ set(SRC_LIBPROJ_CONVERSIONS
   conversions/geoc.cpp
   conversions/geocent.cpp
   conversions/noop.cpp
+  conversions/set.cpp
   conversions/unitconvert.cpp
 )
 
diff --git a/src/pj_list.h b/src/pj_list.h
index 9798a36b..2e1bf759 100644
--- a/src/pj_list.h
+++ b/src/pj_list.h
@@ -134,6 +134,7 @@ PROJ_HEAD(robin, "Robinson")
 PROJ_HEAD(rouss, "Roussilhe Stereographic")
 PROJ_HEAD(rpoly, "Rectangular Polyconic")
 PROJ_HEAD(sch, "Spherical Cross-track Height")
+PROJ_HEAD(set, "Set coordinate value")
 PROJ_HEAD(sinu, "Sinusoidal (Sanson-Flamsteed)")
 PROJ_HEAD(somerc, "Swiss. Obl. Mercator")
 PROJ_HEAD(stere, "Stereographic")