add MeasTable::posArgDeriv and aber*ArgDeriv

Returns derivative of MeasTable::posArg and aber*Arg Polynomial, which is also computed only once and the cached. Avoids unnecessary Polynomial creations in some places which involved several small Array allocations.
casacore · Oct 27, 2015 · fdd1928 · fdd1928
1 parent 00a7605
commit fdd1928
Show file tree

Hide file tree

Showing 4 changed files with 77 additions and 9 deletions.
diff --git a/measures/Measures/Aberration.cc b/measures/Measures/Aberration.cc
@@ -167,9 +167,8 @@ void Aberration::calcAber(Double t) {
     case B1950:
       {
         for (i=0; i<12; i++) {
-          const Polynomial<Double>& aberArgP = MeasTable::aber1950Arg(i);
-          fa(i) = aberArgP(t);
-          dfa(i) = (aberArgP.derivative())(t);
+          fa(i) = MeasTable::aber1950Arg(i)(t);
+          dfa(i) = MeasTable::aber1950ArgDeriv(i)(t);
         }
 	CountedPtr<Matrix<Double> > mul = MeasTable::mulAber1950(t, 1e-6);
         DebugAssert (mul->contiguousStorage(), AipsError);
@@ -213,10 +212,8 @@ void Aberration::calcAber(Double t) {
 	}
       } else {
 	for (i=0; i<13; i++) {
-	  const Polynomial<Double>& aberArgP = MeasTable::aberArg(i);
-
-	  fa(i) = aberArgP(t);
-	  dfa(i) = (aberArgP.derivative())(t);
+	  fa(i) = MeasTable::aberArg(i)(t);
+	  dfa(i) = MeasTable::aberArgDeriv(i)(t);
 	}
 	CountedPtr<Matrix<Double> > mul = MeasTable::mulAber(t, 1e-6);
         DebugAssert (mul->contiguousStorage(), AipsError);

diff --git a/measures/Measures/MeasTable.cc b/measures/Measures/MeasTable.cc
@@ -3236,6 +3236,28 @@ const Polynomial<Double> &MeasTable::aberArg(uInt which) {
   return polyArray[which];
 }
 
+// Derivative aber
+const Polynomial<Double> &MeasTable::aberArgDeriv(uInt which) {
+  static volatile Bool needInit = True;
+  static Polynomial<Double> polyArray[13];
+
+  if (needInit) {
+	const Polynomial<Double> * polyArray_ptrs[13];
+    for (int i=0; i<13; i++) {
+	  polyArray_ptrs[i] = &aberArg(i);
+    }
+    ScopedMutexLock locker(theirMutex);
+    if (needInit) {
+      for (int i=0; i<13; i++) {
+        polyArray[i] = polyArray_ptrs[i]->derivative();
+      }
+      needInit = False;
+    }
+  }
+  DebugAssert(which < 13, AipsError);
+  return polyArray[which];
+}
+
 const Polynomial<Double> &MeasTable::aber1950Arg(uInt which) {
   static volatile Bool needInit = True;
   static Polynomial<Double> polyArray[12];
@@ -3271,6 +3293,28 @@ const Polynomial<Double> &MeasTable::aber1950Arg(uInt which) {
   return polyArray[which];
 }
 
+// Derivative aber1950
+const Polynomial<Double> &MeasTable::aber1950ArgDeriv(uInt which) {
+  static volatile Bool needInit = True;
+  static Polynomial<Double> polyArray[12];
+
+  if (needInit) {
+	const Polynomial<Double> * polyArray_ptrs[12];
+    for (int i=0; i<12; i++) {
+	  polyArray_ptrs[i] = &aber1950Arg(i);
+    }
+    ScopedMutexLock locker(theirMutex);
+    if (needInit) {
+      for (int i=0; i<12; i++) {
+        polyArray[i] = polyArray_ptrs[i]->derivative();
+      }
+      needInit = False;
+    }
+  }
+  DebugAssert(which < 12, AipsError);
+  return polyArray[which];
+}
+
 const Double* MeasTable::mulAberArg(uInt which) {
   static const Double ABERARG[80][6] = {
     {	0,	0,	1,	0,	0,	0},
@@ -3827,6 +3871,29 @@ const Polynomial<Double> &MeasTable::posArg(uInt which) {
                                       POSFUND[i][j]*C::degree);
         }
       }
+
+      needInit = False;
+    }
+  }
+  DebugAssert(which < 12, AipsError);
+  return polyArray[which];
+}
+
+// Derivative of Earth and Sun position polynomial
+const Polynomial<Double> &MeasTable::posArgDeriv(uInt which) {
+  static volatile Bool needInit = True;
+  static Polynomial<Double> polyArray[12];
+
+  if (needInit) {
+	const Polynomial<Double> * polyArray_ptrs[12];
+    for (int i=0; i<12; i++) {
+	  polyArray_ptrs[i] = &posArg(i);
+    }
+    ScopedMutexLock locker(theirMutex);
+    if (needInit) {
+      for (int i=0; i<12; i++) {
+        polyArray[i] = polyArray_ptrs[i]->derivative();
+      }
       needInit = False;
     }
   }

diff --git a/measures/Measures/MeasTable.h b/measures/Measures/MeasTable.h
@@ -289,7 +289,9 @@ class MeasTable {
   // (B1950). 
   // <group>
   static const Polynomial<Double> &aberArg(uInt which);
+  static const Polynomial<Double> &aberArgDeriv(uInt which);
   static const Polynomial<Double> &aber1950Arg(uInt which);
+  static const Polynomial<Double> &aber1950ArgDeriv(uInt which);
   // </group>
 
   // Generate the 'which' vector of the aberration series arguments
@@ -344,6 +346,8 @@ class MeasTable {
   // Fundamental arguments for Soma et al. methods
   // <group>
   static const Polynomial<Double> &posArg(uInt which);
+  // Precomputed derivative of PosArg
+  static const Polynomial<Double> &posArgDeriv(uInt which);
   // </group>
   // Generate the which' vector of the position series arguments
   // <group>

diff --git a/measures/Measures/SolarPos.cc b/measures/Measures/SolarPos.cc
@@ -226,7 +226,7 @@ void SolarPos::calcEarth(Double t) {
 	      } else {
 		for (i=0; i<12; i++) {
 		  fa(i) = MeasTable::posArg(i)(t);
-		  dfa(i) = (MeasTable::posArg(i).derivative())(t);
+		  dfa(i) = MeasTable::posArgDeriv(i)(t);
 		}
                 CountedPtr<Matrix<Double> > mul = MeasTable::mulPosEarthXY(t, 1e-6);
                 DebugAssert (mul->contiguousStorage(), AipsError);
@@ -301,7 +301,7 @@ void SolarPos::calcSun(Double t) {
               } else {
 		for (i=0; i<12; i++) {
 		  fa(i) = MeasTable::posArg(i)(t);
-		  dfa(i) = (MeasTable::posArg(i).derivative())(t);
+		  dfa(i) = MeasTable::posArgDeriv(i)(t);
 		}
                 CountedPtr<Matrix<Double> > mul = MeasTable::mulPosSunXY(t, 1e-6);
                 DebugAssert (mul->contiguousStorage(), AipsError);