Commit 58347e6b authored by Jeannie Backer's avatar Jeannie Backer Committed by Makayla Shepherd
Browse files

Adding SymmetricMatrix and UpperSymmetricMatrix typedefs to LinearAlgebra.

parent e52e88ff
Loading
Loading
Loading
Loading
+12 −1
Original line number Diff line number Diff line
@@ -26,8 +26,9 @@
#include <iostream>

// boost library
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/symmetric.hpp>

// Qt Library
#include <QDebug>
@@ -112,6 +113,16 @@ namespace Isis {
       * API. 
       */
      typedef boost::numeric::ublas::matrix<double> Matrix;
      /**
       * Definition for an Isis::LinearAlgebra::UpperSymmetricMatrix 
       * of doubles. This is a typedef for a boost matrix. 
       */
      typedef boost::numeric::ublas::symmetric_matrix<double> SymmetricMatrix;
      /**
       * Definition for an Isis::LinearAlgebra::UpperSymmetricMatrix 
       * of doubles. This is a typedef for a boost matrix. 
       */
      typedef boost::numeric::ublas::symmetric_matrix<double, upper> UpperSymmetricMatrix;
      /**
       * Definition for an Isis::LinearAlgebra::Vector of doubles. This is a 
       * typedef for a boost vector. 
+21 −21
Original line number Diff line number Diff line
@@ -5,9 +5,9 @@
#include "IException.h"
#include "IString.h"
#include "Latitude.h"
#include "LinearAlgebra.h"
#include "Longitude.h"

using namespace boost::numeric::ublas;
using namespace std;

namespace Isis {
@@ -70,14 +70,14 @@ namespace Isis {
    }

    if(other.p_rectCovar) {
      p_rectCovar = new symmetric_matrix<double, upper>(*other.p_rectCovar);
      p_rectCovar = new LinearAlgebra::UpperSymmetricMatrix(*other.p_rectCovar);
    }
    else {
      p_rectCovar = NULL;
    }

    if(other.p_sphereCovar) {
      p_sphereCovar = new symmetric_matrix<double, upper>(*other.p_sphereCovar);
      p_sphereCovar = new LinearAlgebra::UpperSymmetricMatrix(*other.p_sphereCovar);
    }
    else {
      p_sphereCovar = NULL;
@@ -132,7 +132,7 @@ namespace Isis {
   *
   */
  SurfacePoint::SurfacePoint(const Latitude &lat, const Longitude &lon,
      const Distance &radius, const symmetric_matrix<double, upper> &covar) {
      const Distance &radius, const LinearAlgebra::UpperSymmetricMatrix &covar) {
    InitCovariance();
    InitPoint();
    InitRadii();
@@ -190,7 +190,7 @@ namespace Isis {
   * @param covar  The variance/covariance matrix of the point
   */
  SurfacePoint::SurfacePoint(const Displacement &x, const Displacement &y,
      const Displacement &z, const symmetric_matrix<double, upper> &covar) {
      const Displacement &z, const LinearAlgebra::UpperSymmetricMatrix &covar) {
    InitCovariance();
    InitPoint();
    InitRadii();
@@ -318,7 +318,7 @@ namespace Isis {
   * @return void
   */
  void SurfacePoint::SetRectangular(Displacement x, Displacement y, Displacement z,
                                    const symmetric_matrix<double,upper>& covar) {
                                    const LinearAlgebra::UpperSymmetricMatrix &covar) {
    SetRectangularPoint(x, y, z);
    SetRectangularMatrix(covar);
  }
@@ -342,7 +342,7 @@ namespace Isis {
      throw IException(IException::User, msg, _FILEINFO_);
    }

    symmetric_matrix<double,upper> covar(3);
    LinearAlgebra::UpperSymmetricMatrix covar(3);
    covar.clear();
    covar(0,0) = xSigma.meters() * xSigma.meters();
    covar(1,1) = ySigma.meters() * ySigma.meters();
@@ -359,7 +359,7 @@ namespace Isis {
   * @return void
   */
  void SurfacePoint::SetRectangularMatrix(
       const symmetric_matrix<double, upper> &covar) {
       const LinearAlgebra::UpperSymmetricMatrix &covar) {
    // Make sure the point is set first
    if (!Valid()) {
      IString msg = "A point must be set before a variance/covariance matrix "
@@ -371,7 +371,7 @@ namespace Isis {
      *p_rectCovar = covar;
    }
    else {
      p_rectCovar = new symmetric_matrix<double, upper>(covar);
      p_rectCovar = new LinearAlgebra::UpperSymmetricMatrix(covar);
    }

    SpiceDouble rectMat[3][3];
@@ -406,7 +406,7 @@ namespace Isis {
    J[2][2] = p_z->meters() / radius;

    if(!p_sphereCovar)
      p_sphereCovar = new symmetric_matrix<double, upper>(3);
      p_sphereCovar = new LinearAlgebra::UpperSymmetricMatrix(3);

    SpiceDouble mat[3][3];
    mxm_c (J, rectMat, mat);
@@ -486,7 +486,7 @@ namespace Isis {
   * @param covar Spherical variance/covariance matrix in m*m
   */
  void SurfacePoint::SetSpherical(const Latitude &lat, const Longitude &lon,
      const Distance &radius, const symmetric_matrix<double, upper> &covar) {
      const Distance &radius, const LinearAlgebra::UpperSymmetricMatrix &covar) {
    SetSphericalPoint(lat, lon, radius);
    SetSphericalMatrix(covar);
  }
@@ -518,7 +518,7 @@ namespace Isis {
                                        const Angle &lonSigma,
                                        const Distance &radiusSigma) {
    if (latSigma.isValid() && lonSigma.isValid() && radiusSigma.isValid()) {
      symmetric_matrix<double,upper> covar(3);
      LinearAlgebra::UpperSymmetricMatrix covar(3);
      covar.clear();

      double sphericalCoordinate;
@@ -584,7 +584,7 @@ namespace Isis {
   * @return void
   */
  void SurfacePoint::SetSphericalMatrix(
     const symmetric_matrix<double, upper> & covar) {
     const LinearAlgebra::UpperSymmetricMatrix &covar) {

    // Make sure the point is set first
    if (!Valid()) {
@@ -597,7 +597,7 @@ namespace Isis {
      *p_sphereCovar = covar;
    }
    else {
      p_sphereCovar = new symmetric_matrix<double, upper>(covar);
      p_sphereCovar = new LinearAlgebra::UpperSymmetricMatrix(covar);
    }

    SpiceDouble sphereMat[3][3];
@@ -637,7 +637,7 @@ namespace Isis {
    J[2][2] = sinPhi;

    if(!p_rectCovar)
      p_rectCovar = new symmetric_matrix<double, upper>(3);
      p_rectCovar = new LinearAlgebra::UpperSymmetricMatrix(3);

    SpiceDouble mat[3][3];
    mxm_c (J, sphereMat, mat);
@@ -825,10 +825,10 @@ namespace Isis {
  }


  symmetric_matrix<double, upper> SurfacePoint::GetRectangularMatrix()
  LinearAlgebra::UpperSymmetricMatrix SurfacePoint::GetRectangularMatrix()
      const {
    if(!p_rectCovar) {
      symmetric_matrix<double, upper> tmp(3);
      LinearAlgebra::UpperSymmetricMatrix tmp(3);
      tmp.clear();
      return tmp;
    }
@@ -977,9 +977,9 @@ namespace Isis {
  }


  symmetric_matrix<double, upper> SurfacePoint::GetSphericalMatrix() const {
  LinearAlgebra::UpperSymmetricMatrix SurfacePoint::GetSphericalMatrix() const {
    if(!p_sphereCovar) {
      symmetric_matrix<double, upper> tmp(3);
      LinearAlgebra::UpperSymmetricMatrix tmp(3);
      tmp.clear();
      return tmp;
    }
@@ -1192,11 +1192,11 @@ namespace Isis {
      }

      if(other.p_rectCovar) {
        p_rectCovar = new symmetric_matrix<double, upper>(*other.p_rectCovar);
        p_rectCovar = new LinearAlgebra::UpperSymmetricMatrix(*other.p_rectCovar);
      }

      if(other.p_sphereCovar) {
        p_sphereCovar = new symmetric_matrix<double, upper>(*other.p_sphereCovar);
        p_sphereCovar = new LinearAlgebra::UpperSymmetricMatrix(*other.p_sphereCovar);
      }
    }

+10 −20
Original line number Diff line number Diff line
@@ -25,9 +25,6 @@
#include <vector>
#include <cmath>

#include "boost/numeric/ublas/symmetric.hpp"
#include "boost/numeric/ublas/io.hpp"

#include "Displacement.h"
#include "Distance.h"
#include "Angle.h"
@@ -96,8 +93,7 @@ namespace Isis {
          const Distance &radiusSigma);
      SurfacePoint(const Latitude &lat, const Longitude &lon,
                   const Distance &radius,
                   const boost::numeric::ublas::symmetric_matrix
                     <double,boost::numeric::ublas::upper>& covar);
                   const LinearAlgebra::UpperSymmetricMatrix &covar);
      SurfacePoint(const Displacement &x, const Displacement &y,
                   const Displacement &z);
      SurfacePoint(const Displacement &x, const Displacement &y,
@@ -105,8 +101,7 @@ namespace Isis {
          const Distance &zSigma);
      SurfacePoint(const Displacement &x, const Displacement &y,
                   const Displacement &z,
                   const boost::numeric::ublas::symmetric_matrix
                     <double,boost::numeric::ublas::upper>& covar);
                   const LinearAlgebra::UpperSymmetricMatrix &covar);
      ~SurfacePoint();

// Rectangular loading utilities
@@ -115,7 +110,7 @@ namespace Isis {
          const Distance &ySigma=Distance(), const Distance &zSigma=Distance());

      void SetRectangular(const Displacement x, const Displacement y, const Displacement z,
        const boost::numeric::ublas::symmetric_matrix<double,boost::numeric::ublas::upper>& covar);
        const LinearAlgebra::UpperSymmetricMatrix &covar);

      //! Set surface point and sigmas in rectangular coordinates and convert to planetocentric
      void SetRectangularSigmas(const Distance &xSigma, const Distance &ySigma,
@@ -123,7 +118,7 @@ namespace Isis {


      void SetRectangularMatrix(
        const boost::numeric::ublas::symmetric_matrix<double,boost::numeric::ublas::upper>& covar);
        const LinearAlgebra::UpperSymmetricMatrix &covar);

// Spherical loading utilities

@@ -136,14 +131,13 @@ namespace Isis {

      void SetSpherical (const Latitude &lat, const Longitude &lon,
          const Distance &radius,
          const boost::numeric::ublas::symmetric_matrix
            <double,boost::numeric::ublas::upper>& covar);
          const LinearAlgebra::UpperSymmetricMatrix &covar);

      void SetSphericalCoordinates(const Latitude &lat, const Longitude &lon,
                                   const Distance &radius);

      void SetSphericalMatrix(
        const boost::numeric::ublas::symmetric_matrix<double,boost::numeric::ublas::upper>& covar);
        const LinearAlgebra::UpperSymmetricMatrix &covar);

      void SetSphericalSigmas(const Angle &latSigma, const Angle &lonSigma,
                              const Distance &radiusSigma);
@@ -165,8 +159,7 @@ namespace Isis {
      Distance GetXSigma() const;
      Distance GetYSigma() const;
      Distance GetZSigma() const;
      boost::numeric::ublas::symmetric_matrix<double,boost::numeric::ublas::upper> 
        GetRectangularMatrix() const;
      LinearAlgebra::UpperSymmetricMatrix GetRectangularMatrix() const;
      Latitude GetLatitude() const; 
      Longitude GetLongitude() const;
      Distance GetLocalRadius() const;
@@ -178,8 +171,7 @@ namespace Isis {
      double GetLonWeight() const;
      Distance GetLocalRadiusSigma() const;
      double GetLocalRadiusWeight() const;
      boost::numeric::ublas::symmetric_matrix
          <double,boost::numeric::ublas::upper> GetSphericalMatrix() const;
      LinearAlgebra::UpperSymmetricMatrix GetSphericalMatrix() const;

// Computational methods
      Distance GetDistanceToPoint(const SurfacePoint &other) const;
@@ -209,11 +201,9 @@ namespace Isis {
      Displacement *p_y;
      Displacement *p_z;
      //! 3x3 upper triangular covariance matrix rectangular coordinates
      boost::numeric::ublas::symmetric_matrix
          <double,boost::numeric::ublas::upper> *p_rectCovar;
      LinearAlgebra::UpperSymmetricMatrix *p_rectCovar;
      //! 3x3 upper triangular covariance matrix ocentric coordinates
      boost::numeric::ublas::symmetric_matrix
          <double,boost::numeric::ublas::upper> *p_sphereCovar;
      LinearAlgebra::UpperSymmetricMatrix *p_sphereCovar;
  };
};

+7 −9
Original line number Diff line number Diff line
@@ -8,15 +8,13 @@
#include "Distance.h"
#include "IException.h"
#include "Latitude.h"
#include "LinearAlgebra.h"
#include "Longitude.h"
#include "Preference.h"
#include "SurfacePoint.h"

#include "boost/numeric/ublas/symmetric.hpp"

using namespace Isis;
using namespace std;
using namespace boost::numeric::ublas;

/**
 *
@@ -28,7 +26,7 @@ using namespace boost::numeric::ublas;
 */
int main(int argc, char *argv[]) {
  Isis::Preference::Preferences(true);
  symmetric_matrix<double,upper> cvRect,cvOc;
  LinearAlgebra::UpperSymmetricMatrix cvRect,cvOc;

  try {
    cout << "UnitTest for SurfacePoint" << endl << endl;
@@ -41,7 +39,7 @@ int main(int argc, char *argv[]) {
                   Distance(1000., Distance::Meters),
                   Distance(1000., Distance::Meters));

    symmetric_matrix<double,upper> covar;
    LinearAlgebra::UpperSymmetricMatrix covar;
    covar.resize(3);
    covar.clear();
    covar(0,0) = 100.;
@@ -58,7 +56,7 @@ int main(int argc, char *argv[]) {
    double latSig = spRec.GetLatSigma().degrees();
    double lonSig = spRec.GetLonSigma().degrees();
    double radSig = spRec.GetLocalRadiusSigma().meters();
    symmetric_matrix<double,upper> covarSphere(3);
    LinearAlgebra::UpperSymmetricMatrix covarSphere(3);
    covarSphere.clear();
    covarSphere = spRec.GetSphericalMatrix();
    
@@ -108,7 +106,7 @@ int main(int argc, char *argv[]) {
                          Longitude(lon, Angle::Degrees),
                          Distance(radius, Distance::Meters),
                          covarSphere );
    symmetric_matrix<double,upper> covarRec(3);
    LinearAlgebra::UpperSymmetricMatrix covarRec(3);
    covarRec.clear();
    covarRec = spSphere.GetRectangularMatrix();

@@ -155,7 +153,7 @@ int main(int argc, char *argv[]) {
    double latSig = (spRec.GetLatSigma()).degrees();
    double lonSig = (spRec.GetLonSigma()).degrees();
    double radSig = spRec.GetLocalRadiusSigma().meters();
    symmetric_matrix<double,upper> covarSphere(3);
    LinearAlgebra::UpperSymmetricMatrix covarSphere(3);
    covarSphere.clear();
    covarSphere = spRec.GetSphericalMatrix();
    
@@ -186,7 +184,7 @@ int main(int argc, char *argv[]) {
                         Angle(1.78752107, Angle::Degrees),
                         Distance(38.454887335682053718134171237789,
                                  Distance::Meters));
    symmetric_matrix<double,upper> covarRec(3);
    LinearAlgebra::UpperSymmetricMatrix covarRec(3);
    covarRec.clear();
    covarRec = spSphere1.GetRectangularMatrix();
    spSphere2.SetSpherical(Latitude(0.55850536, Angle::Radians),
+12 −11
Original line number Diff line number Diff line
@@ -37,6 +37,7 @@
#include "ImageList.h"
#include "iTime.h"
#include "Latitude.h"
#include "LinearAlgebra.h"
#include "Longitude.h"
#include "MaximumLikelihoodWFunctions.h"
#include "SpecialPixel.h"
@@ -1030,7 +1031,7 @@ namespace Isis {
    static LinearAlgebra::Matrix coeffImage;
    static LinearAlgebra::Matrix coeffPoint3D(2, 3);
    static LinearAlgebra::Vector coeffRHS(2);
    static boost::numeric::ublas::symmetric_matrix<double, upper> N22(3);
    static LinearAlgebra::UpperSymmetricMatrix N22(3);
    SparseBlockColumnMatrix N12;
    static LinearAlgebra::Vector n2(3);
    boost::numeric::ublas::compressed_vector<double> n1(m_rank);
@@ -1147,7 +1148,7 @@ namespace Isis {
   *
   * @see BundleAdjust::formNormalEquations
   */
  bool BundleAdjust::formMeasureNormals(symmetric_matrix<double, upper>&N22,
  bool BundleAdjust::formMeasureNormals(LinearAlgebra::UpperSymmetricMatrix &N22,
                                        SparseBlockColumnMatrix &N12,
                                        compressed_vector<double> &n1,
                                        vector<double> &n2,
@@ -1157,7 +1158,7 @@ namespace Isis {
                                        vector<double> &coeffRHS,
                                        int observationIndex) {

    static symmetric_matrix<double, upper> N11;
    static LinearAlgebra::UpperSymmetricMatrix N11;
    static matrix<double> N11TargetImage;

    int blockIndex = observationIndex;
@@ -1283,7 +1284,7 @@ namespace Isis {
   *
   * @see BundleAdjust::formNormalEquations
   */
  bool BundleAdjust::formPointNormals(symmetric_matrix<double, upper>&N22,
  bool BundleAdjust::formPointNormals(LinearAlgebra::UpperSymmetricMatrix &N22,
                                      SparseBlockColumnMatrix &N12,
                                      vector<double> &n2,
                                      vector<double> &nj,
@@ -1454,7 +1455,7 @@ namespace Isis {
   *
   * @see BundleAdjust::formPointNormals
   */
  bool BundleAdjust::productATransB(symmetric_matrix <double,upper> &N22,
  bool BundleAdjust::productATransB(LinearAlgebra::UpperSymmetricMatrix &N22,
                                    SparseBlockColumnMatrix &N12,
                                    SparseBlockRowMatrix &Q) {

@@ -1795,11 +1796,11 @@ namespace Isis {
   *
   * @TODO Move to LinearAlgebra
   */
  bool BundleAdjust::invert3x3(symmetric_matrix<double, upper> &m) {
  bool BundleAdjust::invert3x3(LinearAlgebra::UpperSymmetricMatrix &m) {
    double det;
    double den;

    boost::numeric::ublas::symmetric_matrix< double, upper > c = m;
    LinearAlgebra::UpperSymmetricMatrix c = m;

    den = m(0, 0) * (m(1, 1) * m(2, 2) - m(1, 2) * m(2, 1))
          - m(0, 1) * (m(1, 0) * m(2, 2) - m(1, 2) * m(2, 0))
@@ -2657,8 +2658,8 @@ namespace Isis {
    printf("     Time: %s\n\n", currentTime.c_str());

    // create and initialize array of 3x3 matrices for all object points
    std::vector< symmetric_matrix<double> > pointCovariances(numObjectPoints,
                                                             symmetric_matrix<double>(3));
    std::vector< LinearAlgebra::SymmetricMatrix > pointCovariances(numObjectPoints,
                                                                   SymmetricMatrix(3));
    for (int d = 0; d < numObjectPoints; d++) {
      pointCovariances[d].clear();
    }
@@ -2808,7 +2809,7 @@ namespace Isis {
        T.clear();

        // get corresponding point covariance matrix
        boost::numeric::ublas::symmetric_matrix<double> &covariance = pointCovariances[pointIndex];
        LinearAlgebra::SymmetricMatrix &covariance = pointCovariances[pointIndex];

        // get firstQBlock - index i is the key into Q for firstQBlock
        LinearAlgebra::Matrix *firstQBlock = Q.value(i);
@@ -2897,7 +2898,7 @@ namespace Isis {
      }

      // get corresponding point covariance matrix
      boost::numeric::ublas::symmetric_matrix<double> &covariance = pointCovariances[pointIndex];
      LinearAlgebra::SymmetricMatrix &covariance = pointCovariances[pointIndex];

      // Ask Ken what is happening here...Setting just the sigmas is not very accurate
      // Shouldn't we be updating and setting the matrix???  TODO
Loading