From 8eb0f9f2ea15bc3973d1530bc20b33aace864e56 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Marek=20Ne=C4=8Dada?= <marek@necada.org>
Date: Sun, 14 Jul 2019 23:26:03 +0300
Subject: [PATCH] Convertors between csph_t and ccart3_t

Former-commit-id: 0d081dc430d23d8431c8b22d2d9a22de5ba86695
---
 qpms/vectors.h | 38 ++++++++++++++++++++++++++++++++++++++
 1 file changed, 38 insertions(+)

diff --git a/qpms/vectors.h b/qpms/vectors.h
index 9bb94bb..f6abed5 100644
--- a/qpms/vectors.h
+++ b/qpms/vectors.h
@@ -287,6 +287,44 @@ static inline csphvec_t ccart2csphvec(const ccart3_t cartvec, const sph_t at) {
 	return res;
 }
 
+/// Convert sph_t to csph_t.
+static inline csph_t sph2csph(sph_t s) {
+	csph_t cs = {s.r, s.theta, s.phi};
+	return cs;
+}
+
+/// Lossy coordinate transform of ccart3_t to csph_t.
+/** The angle and real part of the radial coordinate are determined
+ *  from the real components of \a \cart. The imaginary part of the radial
+ *  coordinate is then determined as the length of the imaginary 
+ *  part of \a cart *projected onto* the real part of \a cart.
+ *
+ *  N.B. this obviously makes not much sense for purely imaginary vectors
+ *  (and will cause NANs). TODO handle this better.
+ */
+static inline csph_t ccart2csph(const ccart3_t cart) {
+	cart3_t rcart = {creal(cart.x), creal(cart.y), creal(cart.z)};
+	cart3_t icart = {cimag(cart.x), cimag(cart.y), cimag(cart.z)};
+	csph_t sph = sph2csph(cart2sph(rcart));
+	sph.r += I * cart3_dot(icart,rcart) / cart3norm(rcart);
+	return sph;
+}
+
+/// Coordinate transform of csph_t to ccart3_t
+static inline ccart3_t csph2ccart(const csph_t sph) {
+	ccart3_t cart;
+	double sin_th = 
+#ifdef QPMS_VECTORS_NICE_TRANSFORMATIONS
+	   (sph.theta == M_PI) ? 0 :
+#endif
+		sin(sph.theta);
+	cart.x = sph.r * sin_th * cos(sph.phi);
+	cart.y = sph.r * sin_th * sin(sph.phi);
+	cart.z = sph.r * cos(sph.theta);
+	return cart;
+}
+
+
 void print_csphvec(csphvec_t);
 void print_ccart3(ccart3_t);
 void print_cart3(cart3_t);