1 files changed, 22 insertions, 0 deletions
diff --git a/src/pj_qsfn.cpp b/src/pj_qsfn.cpp
new file mode 100644
index 00000000..c18a7b95
--- /dev/null
+++ b/src/pj_qsfn.cpp
@@ -0,0 +1,22 @@
+/* determine small q */
+#include <math.h>
+#include "projects.h"
+
+# define EPSILON 1.0e-7
+
+double pj_qsfn(double sinphi, double e, double one_es) {
+    double con, div1, div2;
+
+    if (e >= EPSILON) {
+        con = e * sinphi;
+        div1 = 1.0 - con * con;
+        div2 = 1.0 + con;
+
+        /* avoid zero division, fail gracefully */
+        if (div1 == 0.0 || div2 == 0.0)
+            return HUGE_VAL;
+
+        return (one_es * (sinphi / div1 - (.5 / e) * log ((1. - con) / div2 )));
+    } else
+        return (sinphi + sinphi);
+}