From b98f386054899a17c12336cb5888394c287172ee Mon Sep 17 00:00:00 2001
From: William Joye <wjoye@cfa.harvard.edu>
Date: Sat, 16 Dec 2017 16:50:32 -0500
Subject: upgrade AST

---
 tksao/frame/fitsimage.C | 101 ++++++++++--------------------------------------
 1 file changed, 21 insertions(+), 80 deletions(-)

diff --git a/tksao/frame/fitsimage.C b/tksao/frame/fitsimage.C
index 090b8d7..cdb2bc6 100644
--- a/tksao/frame/fitsimage.C
+++ b/tksao/frame/fitsimage.C
@@ -3030,67 +3030,22 @@ Coord::Orientation FitsImage::getWCSOrientation(Coord::CoordSystem sys,
   
   setWCSSkyFrame(sys, sky);
 
-  double ang =0;
-  Coord::Orientation rr = Coord::NORMAL;
-  {
-    Vector in[3];
-    Vector out[3];
-    in[0] = center();
-    in[1] = center()+Vector(0,1);
-    in[2] = center()+Vector(1,0);
-    wcsTran(3, in, 1, out);
-    ang = wcsAngle(out[0],out[1],out[2]);
-
-    //    cerr << ang << ':' << radToDeg(ang) << ':';
-    if (!(isnan(ang)||isinf(ang)||(ang == -DBL_MAX)||(ang == DBL_MAX))) {
-      if (hasWCSCel(sys))
-	rr = (ang>=0 ? Coord::NORMAL : Coord::XX);
-      else
-	rr = (ang<=0 ? Coord::NORMAL : Coord::XX);
-    }
-    //    cerr << rr << '=';
-  }
-  if (0) {
-    Vector in[3];
-    Vector out[3];
-    in[0] = center();
-    in[1] = center()+Vector(0,1);
-    in[2] = center()+Vector(1,0);
-    wcsTran(3, in, 1, out);
-    double ang = wcsAngle(out[1],out[0],out[2]);
-
-    //    cerr << ang << ':' << radToDeg(ang) << ':';
-    if (!(isnan(ang)||isinf(ang)||(ang == -DBL_MAX)||(ang == DBL_MAX))) {
-      // backward
-      if (hasWCSCel(sys))
-	rr = (ang<=0 ? Coord::NORMAL : Coord::XX);
-      else
-	rr = (ang>=0 ? Coord::NORMAL : Coord::XX);
-    }
-    cerr << rr << '=';
-  }
-  if (0) {
-    Vector cc = center();
-    Vector wcc = wcsTran(cc, 1);
-    Vector wnorth = wcc + Vector(0,.001);
-    Vector weast = wcc + Vector(.001,0);
-    Vector north = wcsTran(wnorth,0);
-    Vector east = wcsTran(weast,0);
-    Vector nnorth = north.normalize();
-    Vector neast = east.normalize();
-    double w = neast[0]*nnorth[1] - neast[1]*nnorth[0];
-
-    //    cerr << w << ':';
-    if (hasWCSCel(sys))
-      rr = (w<=0 ? Coord::NORMAL : Coord::XX);
-    else
-      rr = (w>=0 ? Coord::NORMAL : Coord::XX);
-
-    cerr << rr;
-  }
-  //  cerr << endl;
-
-  return rr;
+  Vector cc = center();
+  Vector wcc = wcsTran(cc, 1);
+  Vector wnorth = wcc + Vector(0,.001);
+  Vector weast = wcc + Vector(.001,0);
+  Vector north = wcsTran(wnorth,0);
+  Vector east = wcsTran(weast,0);
+  Vector nnorth = north.normalize();
+  Vector neast = east.normalize();
+
+  // take the cross product and see which way the 3rd axis is pointing
+  double ww = neast[0]*nnorth[1] - neast[1]*nnorth[0];
+
+  if (hasWCSCel(sys))
+    return ww<=0 ? Coord::NORMAL : Coord::XX;
+  else
+    return ww>=0 ? Coord::NORMAL : Coord::XX;
 }
 #endif
 
@@ -3126,19 +3081,6 @@ double FitsImage::getWCSRotation(Coord::CoordSystem sys, Coord::SkyFrame sky)
   setWCSSkyFrame(sys, sky);
 
   if (!hasWCSHPX(sys)) {
-    double old =0;
-    {
-      Vector in[2];
-      Vector out[2];
-      in[0] = center();
-      in[1] = center()+Vector(0,1);
-      wcsTran(2, in, 1, out);
-      double ang = wcsAxAngle(out[0], out[1]);
-
-      if (!(isnan(ang)||isinf(ang)||(ang == -DBL_MAX)||(ang == DBL_MAX)))
-	old = getWCSOrientation(sys,sky) == Coord::NORMAL ? ang : -ang;
-    }
-
     Vector in[3];
     Vector out[3];
     in[0] = center();
@@ -3151,11 +3093,9 @@ double FitsImage::getWCSRotation(Coord::CoordSystem sys, Coord::SkyFrame sky)
     if ((!(isnan(ang )||isinf(ang) ||(ang ==-DBL_MAX)||(ang ==DBL_MAX))) &&
 	(!(isnan(ang3)||isinf(ang3)||(ang3==-DBL_MAX)||(ang3==DBL_MAX)))) {
       if ((hasWCSCel(sys) && ang3>0) || (!hasWCSCel(sys) && ang3<0))
-	ang = -ang;
-
-      if (ang != old)
-	cerr << "***BOOM*** " << radToDeg(ang) << '=' << radToDeg(old) << endl;
-      return ang;
+	return -ang;
+      else
+	return ang;
     }
     else
       return 0;
@@ -3173,8 +3113,9 @@ double FitsImage::getWCSRotation(Coord::CoordSystem sys, Coord::SkyFrame sky)
     astSetI(ast_,"Current",current);
     if (!(isnan(ang)||isinf(ang)||(ang == -DBL_MAX)||(ang == DBL_MAX)))
       return ang;
+    else
+      return 0;
   }
-  return 0;
 }
 #endif
 
-- 
cgit v0.12