Isis 3.0 Programmer Reference
Back | Home
BundleObservation.cpp
1 #include "BundleObservation.h"
2 
3 #include <QDebug>
4 #include <QString>
5 #include <QStringList>
6 
7 #include "BundleImage.h"
9 #include "BundleTargetBody.h"
10 #include "Camera.h"
11 #include "LinearAlgebra.h"
12 #include "SpicePosition.h"
13 #include "SpiceRotation.h"
14 
15 using namespace std;
16 
17 namespace Isis {
18 
22  BundleObservation::BundleObservation() {
23  m_serialNumbers.clear();
24  m_imageNames.clear();
25  m_parameterNamesList.clear();
26  m_observationNumber = "";
27  m_instrumentId = "";
28  m_instrumentRotation = NULL;
29  m_instrumentPosition = NULL;
30  m_index = 0;
31  m_weights.clear();
32  m_corrections.clear();
33 // m_solution.clear();
34  m_aprioriSigmas.clear();
35  m_adjustedSigmas.clear();
36  }
37 
38 
48  BundleObservation::BundleObservation(BundleImageQsp image, QString observationNumber,
49  QString instrumentId, BundleTargetBodyQsp bundleTargetBody) {
50  m_serialNumbers.clear();
51  m_imageNames.clear();
52  m_parameterNamesList.clear();
53  m_observationNumber = "";
54  m_instrumentId = "";
55  m_instrumentRotation = NULL;
56  m_instrumentPosition = NULL;
57  m_index = 0;
58  m_weights.clear();
59  m_corrections.clear();
60 // m_solution.clear();
61  m_aprioriSigmas.clear();
62  m_adjustedSigmas.clear();
63 
64  m_observationNumber = observationNumber;
65  m_instrumentId = instrumentId;
66 
67  m_bundleTargetBody = bundleTargetBody;
68 
69  if (image) {
70  append(image);
71  m_serialNumbers.append(image->serialNumber());
72  m_imageNames.append(image->fileName());
73  m_cubeSerialNumberToBundleImageMap.insert(image->serialNumber(), image);
74 
75  // set the observations spice position and rotation objects from the primary image in the
76  // observation (this is, by design at the moment, the first image added to the observation)
77  // if the image, camera, or instrument position/orientation is null, then set to null
78  m_instrumentPosition = (image->camera() ?
79  (image->camera()->instrumentPosition() ?
80  image->camera()->instrumentPosition() : NULL)
81  : NULL);
82  m_instrumentRotation = (image->camera() ?
83  (image->camera()->instrumentRotation() ?
84  image->camera()->instrumentRotation() : NULL)
85  : NULL);
86 
87  // set the observations target body spice rotation object from the primary image in the
88  // observation (this is, by design at the moment, the first image added to the observation)
89  // if the image, camera, or instrument position/orientation is null, then set to null
90 // m_bodyRotation = (image->camera() ?
91 // (image->camera()->bodyRotation() ?
92 // image->camera()->bodyRotation() : NULL)
93 // : NULL);
94  }
95  }
96 
97 
103  BundleObservation::BundleObservation(const BundleObservation &src) {
104  m_serialNumbers = src.m_serialNumbers;
105  m_cubeSerialNumberToBundleImageMap = src.m_cubeSerialNumberToBundleImageMap;
106 
107  m_observationNumber = src.m_observationNumber;
108  m_instrumentId = src.m_instrumentId;
109 
110  m_instrumentPosition = src.m_instrumentPosition;
111  m_instrumentRotation = src.m_instrumentRotation;
112 
113  m_solveSettings = src.m_solveSettings;
114 
115  m_index = src.m_index;
116  }
117 
118 
124  BundleObservation::~BundleObservation() {
125  clear();
126  }
127 
128 
138  BundleObservation &BundleObservation::operator=(const BundleObservation &src) {
139  if (&src != this) {
140  m_serialNumbers = src.m_serialNumbers;
141  m_cubeSerialNumberToBundleImageMap = src.m_cubeSerialNumberToBundleImageMap;
142 
143  m_observationNumber = src.m_observationNumber;
144  m_instrumentId = src.m_instrumentId;
145 
146  m_instrumentPosition = src.m_instrumentPosition;
147  m_instrumentRotation = src.m_instrumentRotation;
148 
149  m_solveSettings = src.m_solveSettings;
150  }
151  return *this;
152  }
153 
154 
164  void BundleObservation::append(const BundleImageQsp &value) {
165  if (value) {
166  m_cubeSerialNumberToBundleImageMap.insert(value->serialNumber(), value);
167  }
169  }
170 
171 
180  BundleImageQsp BundleObservation::imageByCubeSerialNumber(QString cubeSerialNumber) {
181  BundleImageQsp bundleImage;
182 
183  if (m_cubeSerialNumberToBundleImageMap.contains(cubeSerialNumber)) {
184  bundleImage = m_cubeSerialNumberToBundleImageMap.value(cubeSerialNumber);
185  }
186 
187  return bundleImage;
188  }
189 
190 
202  bool BundleObservation::setSolveSettings(BundleObservationSolveSettings solveSettings) {
203  m_solveSettings = BundleObservationSolveSettingsQsp(
204  new BundleObservationSolveSettings(solveSettings));
205 
206  // initialize solution parameters for this observation
207  int nCameraAngleCoefficients = m_solveSettings->numberCameraAngleCoefficientsSolved();
208  int nCameraPositionCoefficients = m_solveSettings->numberCameraPositionCoefficientsSolved();
209 
210  int nParameters = 3*nCameraPositionCoefficients + 2*nCameraAngleCoefficients;
211  if (nCameraAngleCoefficients >= 1 && m_solveSettings->solveTwist()) {
212  nParameters += nCameraAngleCoefficients;
213  }
214  // size vectors and set to zero
215  m_weights.resize(nParameters);
216  m_weights.clear();
217  m_corrections.resize(nParameters);
218  m_corrections.clear();
219 // m_solution.resize(nParameters);
220 // m_solution.clear();
221  m_adjustedSigmas.resize(nParameters);
222  m_adjustedSigmas.clear();
223  m_aprioriSigmas.resize(nParameters);
224  for ( int i = 0; i < nParameters; i++) // initialize apriori sigmas to -1.0
225  m_aprioriSigmas[i] = Isis::Null;
226 
227  if (!initParameterWeights()) {
228  // TODO: some message here!!!!!!!!!!!
229  // TODO: do we need this??? initParameterWeights() never returns false...
230  return false;
231  }
232 
233  return true;
234  }
235 
236 
242  QString BundleObservation::instrumentId() {
243  return m_instrumentId;
244  }
245 
246 
252  SpiceRotation *BundleObservation::spiceRotation() {
253  return m_instrumentRotation;
254  }
255 
256 
262  SpicePosition *BundleObservation::spicePosition() {
263  return m_instrumentPosition;
264  }
265 
266 
272  LinearAlgebra::Vector &BundleObservation::parameterWeights() {
273  return m_weights;
274  }
275 
276 
282  LinearAlgebra::Vector &BundleObservation::parameterCorrections() {
283  return m_corrections;
284  }
285 
286 
291 // LinearAlgebra::Vector &BundleObservation::parameterSolution() {
292 // return m_solution;
293 // }
294 
295 
301  LinearAlgebra::Vector &BundleObservation::aprioriSigmas() {
302  return m_aprioriSigmas;
303  }
304 
305 
311  LinearAlgebra::Vector &BundleObservation::adjustedSigmas() {
312  return m_adjustedSigmas;
313  }
314 
315 
322  const BundleObservationSolveSettingsQsp BundleObservation::solveSettings() {
323  return m_solveSettings;
324  }
325 
326 
335  bool BundleObservation::initializeExteriorOrientation() {
336 
337  if (m_solveSettings->instrumentPositionSolveOption() !=
338  BundleObservationSolveSettings::NoPositionFactors) {
339 
340  double positionBaseTime = 0.0;
341  double positiontimeScale = 0.0;
342  std::vector<double> posPoly1, posPoly2, posPoly3;
343 
344  for (int i = 0; i < size(); i++) {
345  BundleImageQsp image = at(i);
346  SpicePosition *spicePosition = image->camera()->instrumentPosition();
347 
348  if (i > 0) {
349  spicePosition->SetPolynomialDegree(m_solveSettings->spkSolveDegree());
350  spicePosition->SetOverrideBaseTime(positionBaseTime, positiontimeScale);
351  spicePosition->SetPolynomial(posPoly1, posPoly2, posPoly3,
352  m_solveSettings->positionInterpolationType());
353  }
354  else {
355  // first, set the degree of the spk polynomial to be fit for a priori values
356  spicePosition->SetPolynomialDegree(m_solveSettings->spkDegree());
357 
358  // now, set what kind of interpolation to use (SPICE, memcache, hermitecache, polynomial
359  // function, or polynomial function over constant hermite spline)
360  // TODO: verify - I think this actually performs the a priori fit
361  spicePosition->SetPolynomial(m_solveSettings->positionInterpolationType());
362 
363  // finally, set the degree of the spk polynomial actually used in the bundle adjustment
364  spicePosition->SetPolynomialDegree(m_solveSettings->spkSolveDegree());
365 
366  if (m_instrumentPosition) { // ??? TODO: why is this different from rotation code below???
367  positionBaseTime = m_instrumentPosition->GetBaseTime();
368  positiontimeScale = m_instrumentPosition->GetTimeScale();
369  m_instrumentPosition->GetPolynomial(posPoly1, posPoly2, posPoly3);
370  }
371  }
372  }
373  }
374 
375  if (m_solveSettings->instrumentPointingSolveOption() !=
376  BundleObservationSolveSettings::NoPointingFactors) {
377 
378  double rotationBaseTime = 0.0;
379  double rotationtimeScale = 0.0;
380  std::vector<double> anglePoly1, anglePoly2, anglePoly3;
381 
382  for (int i = 0; i < size(); i++) {
383  BundleImageQsp image = at(i);
384  SpiceRotation *spicerotation = image->camera()->instrumentRotation();
385 
386  if (i > 0) {
387  spicerotation->SetPolynomialDegree(m_solveSettings->ckSolveDegree());
388  spicerotation->SetOverrideBaseTime(rotationBaseTime, rotationtimeScale);
389  spicerotation->SetPolynomial(anglePoly1, anglePoly2, anglePoly3,
390  m_solveSettings->pointingInterpolationType());
391  }
392  else {
393  // first, set the degree of the spk polynomial to be fit for a priori values
394  spicerotation->SetPolynomialDegree(m_solveSettings->ckDegree());
395 
396  // now, set what kind of interpolation to use (SPICE, memcache, hermitecache, polynomial
397  // function, or polynomial function over constant hermite spline)
398  // TODO: verify - I think this actually performs the a priori fit
399  spicerotation->SetPolynomial(m_solveSettings->pointingInterpolationType());
400 
401  // finally, set the degree of the spk polynomial actually used in the bundle adjustment
402  spicerotation->SetPolynomialDegree(m_solveSettings->ckSolveDegree());
403 
404  rotationBaseTime = spicerotation->GetBaseTime();
405  rotationtimeScale = spicerotation->GetTimeScale();
406  spicerotation->GetPolynomial(anglePoly1, anglePoly2, anglePoly3);
407  }
408  }
409  }
410 
411  return true;
412  }
413 
414 
420  void BundleObservation::initializeBodyRotation() {
421  std::vector<Angle> raCoefs = m_bundleTargetBody->poleRaCoefs();
422  std::vector<Angle> decCoefs = m_bundleTargetBody->poleDecCoefs();
423  std::vector<Angle> pmCoefs = m_bundleTargetBody->pmCoefs();
424 
425  for (int i = 0; i < size(); i++) {
426  BundleImageQsp image = at(i);
427  image->camera()->bodyRotation()->setPckPolynomial(raCoefs, decCoefs, pmCoefs);
428  }
429  }
430 
431 
438  void BundleObservation::updateBodyRotation() {
439  std::vector<Angle> raCoefs = m_bundleTargetBody->poleRaCoefs();
440  std::vector<Angle> decCoefs = m_bundleTargetBody->poleDecCoefs();
441  std::vector<Angle> pmCoefs = m_bundleTargetBody->pmCoefs();
442 
443  for (int i = 0; i < size(); i++) {
444  BundleImageQsp image = at(i);
445  image->camera()->bodyRotation()->setPckPolynomial(raCoefs, decCoefs, pmCoefs);
446  }
447  }
448 
449 
450 /*
451  bool BundleObservation::initializeExteriorOrientation() {
452 
453  if (m_solveSettings->instrumentPositionSolveOption() !=
454  BundleObservationSolveSettings::NoPositionFactors) {
455 
456  for (int i = 0; i < size(); i++) {
457  BundleImageQsp image = at(i);
458  SpicePosition *spiceposition = image->camera()->instrumentPosition();
459 
460  // first, set the degree of the spk polynomial to be fit for a priori values
461  spiceposition->SetPolynomialDegree(m_solveSettings->spkDegree());
462 
463  // now, set what kind of interpolation to use (SPICE, memcache, hermitecache, polynomial
464  // function, or polynomial function over constant hermite spline)
465  // TODO: verify - I think this actually performs the a priori fit
466  spiceposition->SetPolynomial(m_solveSettings->positionInterpolationType());
467 
468  // finally, set the degree of the spk polynomial actually used in the bundle adjustment
469  spiceposition->SetPolynomialDegree(m_solveSettings->spkSolveDegree());
470  }
471  }
472 
473  if (m_solveSettings->instrumentPointingSolveOption() !=
474  BundleObservationSolveSettings::NoPointingFactors) {
475 
476  for (int i = 0; i < size(); i++) {
477  BundleImageQsp image = at(i);
478  SpiceRotation *spicerotation = image->camera()->instrumentRotation();
479 
480  // first, set the degree of the spk polynomial to be fit for a priori values
481  spicerotation->SetPolynomialDegree(m_solveSettings->ckDegree());
482 
483  // now, set what kind of interpolation to use (SPICE, memcache, hermitecache, polynomial
484  // function, or polynomial function over constant hermite spline)
485  // TODO: verify - I think this actually performs the a priori fit
486  spicerotation->SetPolynomial(m_solveSettings->pointingInterpolationType());
487 
488  // finally, set the degree of the spk polynomial actually used in the bundle adjustment
489  spicerotation->SetPolynomialDegree(m_solveSettings->ckSolveDegree());
490  }
491  }
492 
493  return true;
494  }
495 */
496 
497 
512  bool BundleObservation::initParameterWeights() {
513 
514  // weights for
515  double posWeight = 0.0; // position
516  double velWeight = 0.0; // velocity
517  double accWeight = 0.0; // acceleration
518  double angWeight = 0.0; // angles
519  double angVelWeight = 0.0; // angular velocity
520  double angAccWeight = 0.0; // angular acceleration
521 
522  QList<double> aprioriPointingSigmas = m_solveSettings->aprioriPointingSigmas();
523  QList<double> aprioriPositionSigmas = m_solveSettings->aprioriPositionSigmas();
524 
525  int nCamPosCoeffsSolved = 3 *m_solveSettings->numberCameraPositionCoefficientsSolved();
526 
527  int nCamAngleCoeffsSolved;
528  if (m_solveSettings->solveTwist()) {
529  nCamAngleCoeffsSolved = 3 *m_solveSettings->numberCameraAngleCoefficientsSolved();
530  }
531  else {
532  nCamAngleCoeffsSolved = 2 *m_solveSettings->numberCameraAngleCoefficientsSolved();
533  }
534 
535  if (aprioriPositionSigmas.size() >= 1 && aprioriPositionSigmas.at(0) > 0.0) {
536  posWeight = aprioriPositionSigmas.at(0);
537  posWeight = 1.0 / (posWeight *posWeight * 1.0e-6);
538  }
539  if (aprioriPositionSigmas.size() >= 2 && aprioriPositionSigmas.at(1) > 0.0) {
540  velWeight = aprioriPositionSigmas.at(1);
541  velWeight = 1.0 / (velWeight *velWeight * 1.0e-6);
542  }
543  if (aprioriPositionSigmas.size() >= 3 && aprioriPositionSigmas.at(2) > 0.0) {
544  accWeight = aprioriPositionSigmas.at(2);
545  accWeight = 1.0 / (accWeight *accWeight * 1.0e-6);
546  }
547 
548  if (aprioriPointingSigmas.size() >= 1 && aprioriPointingSigmas.at(0) > 0.0) {
549  angWeight = aprioriPointingSigmas.at(0);
550  angWeight = 1.0 / (angWeight *angWeight * DEG2RAD * DEG2RAD);
551  }
552  if (aprioriPointingSigmas.size() >= 2 && aprioriPointingSigmas.at(1) > 0.0) {
553  angVelWeight = aprioriPointingSigmas.at(1);
554  angVelWeight = 1.0 / (angVelWeight * angVelWeight * DEG2RAD * DEG2RAD);
555  }
556  if (aprioriPointingSigmas.size() >= 3 && aprioriPointingSigmas.at(2) > 0.0) {
557  angAccWeight = aprioriPointingSigmas.at(2);
558  angAccWeight = 1.0 / (angAccWeight * angAccWeight * DEG2RAD * DEG2RAD);
559  }
560 
561  int nSpkTerms = m_solveSettings->spkSolveDegree()+1;
562  nSpkTerms = m_solveSettings->numberCameraPositionCoefficientsSolved();
563  for ( int i = 0; i < nCamPosCoeffsSolved; i++) {
564  if (i % nSpkTerms == 0) {
565  m_aprioriSigmas[i] = aprioriPositionSigmas.at(0);
566  m_weights[i] = posWeight;
567  }
568  if (i % nSpkTerms == 1) {
569  m_aprioriSigmas[i] = aprioriPositionSigmas.at(1);
570  m_weights[i] = velWeight;
571  }
572  if (i % nSpkTerms == 2) {
573  m_aprioriSigmas[i] = aprioriPositionSigmas.at(2);
574  m_weights[i] = accWeight;
575  }
576  }
577 
578  int nCkTerms = m_solveSettings->ckSolveDegree()+1;
579  nCkTerms = m_solveSettings->numberCameraAngleCoefficientsSolved();
580  for ( int i = 0; i < nCamAngleCoeffsSolved; i++) {
581  if (i % nCkTerms == 0) {
582  m_aprioriSigmas[nCamPosCoeffsSolved + i] = aprioriPointingSigmas.at(0);
583  m_weights[nCamPosCoeffsSolved + i] = angWeight;
584  }
585  if (i % nCkTerms == 1) {
586  m_aprioriSigmas[nCamPosCoeffsSolved + i] = aprioriPointingSigmas.at(1);
587  m_weights[nCamPosCoeffsSolved + i] = angVelWeight;
588  }
589  if (i % nCkTerms == 2) {
590  m_aprioriSigmas[nCamPosCoeffsSolved + i] = aprioriPointingSigmas.at(2);
591  m_weights[nCamPosCoeffsSolved + i] = angAccWeight;
592  }
593  }
594 
595 // for ( int i = 0; i < (int)m_weights.size(); i++ )
596 // std::cout << m_weights[i] << std::endl;
597 
598  return true;
599  }
600 
601 
618  bool BundleObservation::applyParameterCorrections(LinearAlgebra::Vector corrections) {
619 
620  int index = 0;
621 
622  try {
623  int nCameraAngleCoefficients = m_solveSettings->numberCameraAngleCoefficientsSolved();
624  int nCameraPositionCoefficients = m_solveSettings->numberCameraPositionCoefficientsSolved();
625 
627  = m_solveSettings->instrumentPositionSolveOption();
628  if (positionOption != BundleObservationSolveSettings::NoPositionFactors) {
629 
630  if (!m_instrumentPosition) {
631  QString msg = "Instrument position is NULL, but position solve option is ";
632  msg.append(BundleObservationSolveSettings::instrumentPositionSolveOptionToString(
633  positionOption));
634  throw IException(IException::Unknown, msg, _FILEINFO_);
635  }
636 
637  std::vector<double> coefX(nCameraPositionCoefficients);
638  std::vector<double> coefY(nCameraPositionCoefficients);
639  std::vector<double> coefZ(nCameraPositionCoefficients);
640 
641  m_instrumentPosition->GetPolynomial(coefX, coefY, coefZ);
642 
643  // update X coordinate coefficient(s) and sum parameter correction
644  for (int i = 0; i < nCameraPositionCoefficients; i++) {
645  coefX[i] += corrections(index);
646  index++;
647  }
648 
649  // update Y coordinate coefficient(s) and sum parameter correction
650  for (int i = 0; i < nCameraPositionCoefficients; i++) {
651  coefY[i] += corrections(index);
652  index++;
653  }
654 
655  // update Z coordinate coefficient(s) and sum parameter correction
656  for (int i = 0; i < nCameraPositionCoefficients; i++) {
657  coefZ[i] += corrections(index);
658  index++;
659  }
660 
661  // apply updates to all images in observation
662  for (int i = 0; i < size(); i++) {
663  BundleImageQsp image = at(i);
664  SpicePosition *spiceposition = image->camera()->instrumentPosition();
665  spiceposition->SetPolynomial(coefX, coefY, coefZ,
666  m_solveSettings->positionInterpolationType());
667  }
668  }
669 
671  = m_solveSettings->instrumentPointingSolveOption();
672  if (pointingOption != BundleObservationSolveSettings::NoPointingFactors) {
673 
674  if (!m_instrumentRotation) {
675  QString msg = "Instrument rotation is NULL, but pointing solve option is ";
676  msg.append(BundleObservationSolveSettings::instrumentPointingSolveOptionToString(
677  pointingOption));
678  throw IException(IException::Unknown, msg, _FILEINFO_);
679  }
680 
681  std::vector<double> coefRA(nCameraPositionCoefficients);
682  std::vector<double> coefDEC(nCameraPositionCoefficients);
683  std::vector<double> coefTWI(nCameraPositionCoefficients);
684 
685  m_instrumentRotation->GetPolynomial(coefRA, coefDEC, coefTWI);
686 
687  // update RA coefficient(s)
688  for (int i = 0; i < nCameraAngleCoefficients; i++) {
689  coefRA[i] += corrections(index);
690  index++;
691  }
692 
693  // update DEC coefficient(s)
694  for (int i = 0; i < nCameraAngleCoefficients; i++) {
695  coefDEC[i] += corrections(index);
696  index++;
697  }
698 
699  if (m_solveSettings->solveTwist()) {
700  // update TWIST coefficient(s)
701  for (int i = 0; i < nCameraAngleCoefficients; i++) {
702  coefTWI[i] += corrections(index);
703  index++;
704  }
705  }
706 
707  // apply updates to all images in observation
708  for (int i = 0; i < size(); i++) {
709  BundleImageQsp image = at(i);
710  SpiceRotation *spiceRotation = image->camera()->instrumentRotation();
711  spiceRotation->SetPolynomial(coefRA, coefDEC, coefTWI,
712  m_solveSettings->pointingInterpolationType());
713  }
714  }
715 
716  // update corrections
717  m_corrections += corrections;
718 
719  }
720  catch (IException &e) {
721  QString msg = "Unable to apply parameter corrections to BundleObservation.";
722  throw IException(e, IException::Unknown, msg, _FILEINFO_);
723  }
724  return true;
725  }
726 
727 
733  int BundleObservation::numberPositionParameters() {
734  return 3.0 * m_solveSettings->numberCameraPositionCoefficientsSolved();
735  }
736 
737 
743  int BundleObservation::numberPointingParameters() {
744  int angleCoefficients = m_solveSettings->numberCameraAngleCoefficientsSolved();
745 
746  if (m_solveSettings->solveTwist()) {
747  return 3.0 * angleCoefficients;
748  }
749  return 2.0 * angleCoefficients;
750  }
751 
752 
761  int BundleObservation::numberParameters() {
762  return numberPositionParameters() + numberPointingParameters();
763  }
764 
765 
771  void BundleObservation::setIndex(int n) {
772  m_index = n;
773  }
774 
775 
781  int BundleObservation::index() {
782  return m_index;
783  }
784 
785 
802  QString BundleObservation::formatBundleOutputString(bool errorPropagation, bool imageCSV) {
803 
804  std::vector<double> coefX;
805  std::vector<double> coefY;
806  std::vector<double> coefZ;
807  std::vector<double> coefRA;
808  std::vector<double> coefDEC;
809  std::vector<double> coefTWI;
810 
811  int nPositionCoefficients = m_solveSettings->numberCameraPositionCoefficientsSolved();
812  int nPointingCoefficients = m_solveSettings->numberCameraAngleCoefficientsSolved();
813 
814  // Indicate if we need to obtain default position or pointing values
815  bool useDefaultPosition = false;
816  bool useDefaultPointing = false;
817  // Indicate if we need to use default values when not solving twist
818  bool useDefaultTwist = !(m_solveSettings->solveTwist());
819 
820  // If we aren't solving for position, set the number of coefficients to 1 so we can output the
821  // instrumentPosition's center coordinate values for X, Y, and Z
822  if (nPositionCoefficients == 0) {
823  nPositionCoefficients = 1;
824  useDefaultPosition = true;
825  }
826  // If we arent' solving for pointing, set the number of coefficients to 1 so we can output the
827  // instrumentPointing's center angles for RA, DEC, and TWI
828  if (nPointingCoefficients == 0) {
829  nPointingCoefficients = 1;
830  useDefaultPointing = true;
831  }
832 
833  // Force number of position and pointing parameters to each be 3 (X,Y,Z; RA,DEC,TWI)
834  // so we can always output a value for them
835  int nPositionParameters = 3 * nPositionCoefficients;
836  int nPointingParameters = 3 * nPointingCoefficients;
837  int nParameters = nPositionParameters + nPointingParameters;
838 
839  coefX.resize(nPositionCoefficients);
840  coefY.resize(nPositionCoefficients);
841  coefZ.resize(nPositionCoefficients);
842  coefRA.resize(nPointingCoefficients);
843  coefDEC.resize(nPointingCoefficients);
844  coefTWI.resize(nPointingCoefficients);
845 
846  if (m_instrumentPosition) {
847  if (!useDefaultPosition) {
848  m_instrumentPosition->GetPolynomial(coefX, coefY, coefZ);
849  }
850  // Use the position's center coordinate if not solving for spacecraft position
851  else {
852  const std::vector<double> centerCoord = m_instrumentPosition->GetCenterCoordinate();
853  coefX[0] = centerCoord[0];
854  coefY[0] = centerCoord[1];
855  coefZ[0] = centerCoord[2];
856  }
857  }
858 
859  if (m_instrumentRotation) {
860  if (!useDefaultPointing) {
861  m_instrumentRotation->GetPolynomial(coefRA, coefDEC, coefTWI);
862  }
863  // Use the pointing's center angles if not solving for pointing (rotation)
864  else {
865  const std::vector<double> centerAngles = m_instrumentRotation->GetCenterAngles();
866  coefRA[0] = centerAngles[0];
867  coefDEC[0] = centerAngles[1];
868  coefTWI[0] = centerAngles[2];
869  }
870  }
871 
872  // for convenience, create vectors of parameters names and values in the correct sequence
873  std::vector<double> finalParameterValues;
874  QStringList parameterNamesList;
875 
876  if (!imageCSV) {
877 
878  QString str("%1(t%2)");
879 
880  if (nPositionCoefficients > 0) {
881  for (int i = 0; i < nPositionCoefficients; i++) {
882  finalParameterValues.push_back(coefX[i]);
883  if (i == 0)
884  parameterNamesList.append( str.arg(" X ").arg("0") );
885  else
886  parameterNamesList.append( str.arg(" ").arg(i) );
887  }
888  for (int i = 0; i < nPositionCoefficients; i++) {
889  finalParameterValues.push_back(coefY[i]);
890  if (i == 0)
891  parameterNamesList.append( str.arg(" Y ").arg("0") );
892  else
893  parameterNamesList.append( str.arg(" ").arg(i) );
894  }
895  for (int i = 0; i < nPositionCoefficients; i++) {
896  finalParameterValues.push_back(coefZ[i]);
897  if (i == 0)
898  parameterNamesList.append( str.arg(" Z ").arg("0") );
899  else
900  parameterNamesList.append( str.arg(" ").arg(i) );
901  }
902  }
903  if (nPointingCoefficients > 0) {
904  for (int i = 0; i < nPointingCoefficients; i++) {
905  finalParameterValues.push_back(coefRA[i] * RAD2DEG);
906  if (i == 0)
907  parameterNamesList.append( str.arg(" RA ").arg("0") );
908  else
909  parameterNamesList.append( str.arg(" ").arg(i) );
910  }
911  for (int i = 0; i < nPointingCoefficients; i++) {
912  finalParameterValues.push_back(coefDEC[i] * RAD2DEG);
913  if (i == 0)
914  parameterNamesList.append( str.arg("DEC ").arg("0") );
915  else
916  parameterNamesList.append( str.arg(" ").arg(i) );
917  }
918  for (int i = 0; i < nPointingCoefficients; i++) {
919  finalParameterValues.push_back(coefTWI[i] * RAD2DEG);
920  if (i == 0)
921  parameterNamesList.append( str.arg("TWI ").arg("0") );
922  else
923  parameterNamesList.append( str.arg(" ").arg(i) );
924  }
925  }
926 
927  }// end if(!imageCSV)
928 
929  else {
930  if (nPositionCoefficients > 0) {
931  for (int i = 0; i < nPositionCoefficients; i++) {
932  finalParameterValues.push_back(coefX[i]);
933  }
934  for (int i = 0; i < nPositionCoefficients; i++) {
935  finalParameterValues.push_back(coefY[i]);
936  }
937  for (int i = 0; i < nPositionCoefficients; i++) {
938  finalParameterValues.push_back(coefZ[i]);
939  }
940  }
941  if (nPointingCoefficients > 0) {
942  for (int i = 0; i < nPointingCoefficients; i++) {
943  finalParameterValues.push_back(coefRA[i] * RAD2DEG);
944  }
945  for (int i = 0; i < nPointingCoefficients; i++) {
946  finalParameterValues.push_back(coefDEC[i] * RAD2DEG);
947  }
948  for (int i = 0; i < nPointingCoefficients; i++) {
949  finalParameterValues.push_back(coefTWI[i] * RAD2DEG);
950  }
951  }
952  }//end else
953 
954  // Save the list of parameter names we've accumulated above
955  m_parameterNamesList = parameterNamesList;
956 
957  QString finalqStr = "";
958  QString qStr = "";
959 
960  // Set up default values when we are using default position
961  QString sigma = "N/A";
962  QString adjustedSigma = "N/A";
963  double correction = 0.0;
964 
965  // this implies we're writing to bundleout.txt
966  if (!imageCSV) {
967  // position parameters
968  for (int i = 0; i < nPositionParameters; i++) {
969  // If not using the default position, we can correctly access sigmas and corrections
970  // members
971  if (!useDefaultPosition) {
972  correction = m_corrections(i);
973  adjustedSigma = QString::number(m_adjustedSigmas[i], 'f', 8);
974  sigma = ( IsSpecial(m_aprioriSigmas[i]) ? "FREE" : toString(m_aprioriSigmas[i], 8) );
975  }
976  if (errorPropagation) {
977  qStr = QString("%1%2%3%4%5%6\n").
978  arg( parameterNamesList.at(i) ).
979  arg(finalParameterValues[i] - correction, 17, 'f', 8).
980  arg(correction, 21, 'f', 8).
981  arg(finalParameterValues[i], 20, 'f', 8).
982  arg(sigma, 18).
983  arg(adjustedSigma, 18);
984  }
985  else {
986  qStr = QString("%1%2%3%4%5%6\n").
987  arg( parameterNamesList.at(i) ).
988  arg(finalParameterValues[i] - correction, 17, 'f', 8).
989  arg(correction, 21, 'f', 8).
990  arg(finalParameterValues[i], 20, 'f', 8).
991  arg(sigma, 18).
992  arg("N/A", 18);
993  }
994  finalqStr += qStr;
995  }
996 
997  // We need to use an offset of -3 (1 coef; X,Y,Z) if we used the default center coordinate
998  // (i.e. we did not solve for position), as m_corrections and m_*sigmas are populated
999  // according to which parameters are solved
1000  int offset = 0;
1001  if (useDefaultPosition) {
1002  offset = 3;
1003  }
1004  // pointing parameters
1005  for (int i = nPositionParameters; i < nParameters; i++) {
1006  if (!useDefaultPointing) {
1007  // If solving camera and not solving for twist, provide default values for twist to
1008  // prevent bad indexing into m_corrections and m_*sigmas
1009  // TWIST is last parameter, which corresponds to nParameters - nPointingCoefficients
1010  if ( (i >= nParameters - nPointingCoefficients) && useDefaultTwist) {
1011  correction = 0.0;
1012  adjustedSigma = "N/A";
1013  sigma = "N/A";
1014  }
1015  else {
1016  correction = m_corrections(i - offset);
1017  adjustedSigma = QString::number(m_adjustedSigmas(i-offset) * RAD2DEG, 'f', 8);
1018  sigma = ( IsSpecial(m_aprioriSigmas[i - offset]) ? "FREE" :
1019  toString(m_aprioriSigmas[i-offset], 8) );
1020  }
1021  }
1022  // We are using default pointing, so provide default correction and sigma values to output
1023  else {
1024  correction = 0.0;
1025  adjustedSigma = "N/A";
1026  sigma = "N/A";
1027  }
1028  if (errorPropagation) {
1029  qStr = QString("%1%2%3%4%5%6\n").
1030  arg( parameterNamesList.at(i) ).
1031  arg( (finalParameterValues[i] - correction * RAD2DEG), 17, 'f', 8).
1032  arg(correction * RAD2DEG, 21, 'f', 8).
1033  arg(finalParameterValues[i], 20, 'f', 8).
1034  arg(sigma, 18).
1035  arg(adjustedSigma, 18);
1036  }
1037  else {
1038  qStr = QString("%1%2%3%4%5%6\n").
1039  arg( parameterNamesList.at(i) ).
1040  arg( (finalParameterValues[i] - correction * RAD2DEG), 17, 'f', 8).
1041  arg(correction * RAD2DEG, 21, 'f', 8).
1042  arg(finalParameterValues[i], 20, 'f', 8).
1043  arg(sigma, 18).
1044  arg("N/A", 18);
1045  }
1046  finalqStr += qStr;
1047  }
1048 
1049  }
1050  // this implies we're writing to images.csv
1051  else {
1052  // position parameters
1053  for (int i = 0; i < nPositionParameters; i++) {
1054  if (!useDefaultPosition) {
1055  correction = m_corrections(i);
1056  adjustedSigma = QString::number(m_adjustedSigmas[i], 'f', 8);
1057  sigma = ( IsSpecial(m_aprioriSigmas[i]) ? "FREE" : toString(m_aprioriSigmas[i], 8) );
1058  }
1059  // Provide default values for position if not solving position
1060  else {
1061  correction = 0.0;
1062  adjustedSigma = "N/A";
1063  sigma = "N/A";
1064  }
1065  qStr = "";
1066  if (errorPropagation) {
1067  qStr += toString(finalParameterValues[i] - correction) + ",";
1068  qStr += toString(correction) + ",";
1069  qStr += toString(finalParameterValues[i]) + ",";
1070  qStr += sigma + ",";
1071  qStr += adjustedSigma + ",";
1072  }
1073  else {
1074  qStr += toString(finalParameterValues[i] - correction) + ",";
1075  qStr += toString(correction) + ",";
1076  qStr += toString(finalParameterValues[i]) + ",";
1077  qStr += sigma + ",";
1078  qStr += "N/A,";
1079  }
1080  finalqStr += qStr;
1081  }
1082 
1083  // If not solving position, we need to offset access to correction and sigma members by -3
1084  // (X,Y,Z) since m_corrections and m_*sigmas are populated according to which parameters are
1085  // solved
1086  int offset = 0;
1087  if (useDefaultPosition) {
1088  offset = 3;
1089  }
1090  // pointing parameters
1091  for (int i = nPositionParameters; i < nParameters; i++) {
1092  if (!useDefaultPointing) {
1093  // Use default values if solving camera but not solving for TWIST to prevent bad indexing
1094  // into m_corrections and m_*sigmas
1095  if ( (i >= nParameters - nPointingCoefficients) && useDefaultTwist) {
1096  correction = 0.0;
1097  adjustedSigma = "N/A";
1098  sigma = "N/A";
1099  }
1100  else {
1101  correction = m_corrections(i - offset);
1102  adjustedSigma = QString::number(m_adjustedSigmas(i-offset) * RAD2DEG, 'f', 8);
1103  sigma = ( IsSpecial(m_aprioriSigmas[i-offset]) ? "FREE" :
1104  toString(m_aprioriSigmas[i-offset], 8) );
1105  }
1106  }
1107  // Provide default values for pointing if not solving pointing
1108  else {
1109  correction = 0.0;
1110  adjustedSigma = "N/A";
1111  sigma = "N/A";
1112  }
1113  qStr = "";
1114  if (errorPropagation) {
1115  qStr += toString(finalParameterValues[i] - correction * RAD2DEG) + ",";
1116  qStr += toString(correction * RAD2DEG) + ",";
1117  qStr += toString(finalParameterValues[i]) + ",";
1118  qStr += sigma + ",";
1119  qStr += adjustedSigma + ",";
1120  }
1121  else {
1122  qStr += toString(finalParameterValues[i] - correction * RAD2DEG) + ",";
1123  qStr += toString(correction * RAD2DEG) + ",";
1124  qStr += toString(finalParameterValues[i]) + ",";
1125  qStr += sigma + ",";
1126  qStr += "N/A,";
1127  }
1128  finalqStr += qStr;
1129  }
1130  }
1131 
1132  return finalqStr;
1133  }
1134 
1135 
1141  QStringList BundleObservation::parameterList() {
1142  return m_parameterNamesList;
1143  }
1144 
1145 
1151  QStringList BundleObservation::imageNames() {
1152  return m_imageNames;
1153  }
1154 }
void SetOverrideBaseTime(double baseTime, double timeScale)
Set an override base time to be used with observations on scanners to allow all images in an observat...
const double Null
Value for an Isis Null pixel.
Definition: SpecialPixel.h:109
QString m_observationNumber
This is typically equivalent to serial number except in the case of &quot;observation mode&quot; (e...
void SetPolynomialDegree(int degree)
Set the degree of the polynomials to be fit to the three camera angles for the time period covered by...
void GetPolynomial(std::vector< double > &abcAng1, std::vector< double > &abcAng2, std::vector< double > &abcAng3)
Return the coefficients of a polynomial fit to each of the three camera angles for the time period co...
QSharedPointer< BundleObservationSolveSettings > BundleObservationSolveSettingsQsp
Definition for BundleObservationSolveSettingsQsp, a QSharedPointer to a &lt; BundleObservationSolveSet...
QMap< QString, BundleImageQsp > m_cubeSerialNumberToBundleImageMap
Map between cube serial number and BundleImage pointers.
const double DEG2RAD(0.017453292519943295769237)
Multiplier for converting from degrees to radians.
QSharedPointer< BundleTargetBody > BundleTargetBodyQsp
Definition for BundleTargetBodyQsp, a QSharedPointer to a BundleTargetBody.
Class for bundle observations.
QString toString(bool boolToConvert)
Global function to convert a boolean to a string.
Definition: IString.cpp:226
BundleObservationSolveSettingsQsp m_solveSettings
Solve settings for this observation.
boost::numeric::ublas::vector< double > Vector
Definition for an Isis::LinearAlgebra::Vector of doubles.
double GetBaseTime()
Accessor method to get the rotation base time.
SpiceRotation * m_instrumentRotation
Instrument spice rotation (in primary image).
#define _FILEINFO_
Macro for the filename and line number.
Definition: IException.h:38
void SetPolynomialDegree(int degree)
Set the polynomial degree.
void SetOverrideBaseTime(double baseTime, double timeScale)
Set an override base time to be used with observations on scanners to allow all images in an observat...
QStringList m_serialNumbers
List of all cube serial numbers in observation.
bool IsSpecial(const double d)
Returns if the input pixel is special.
Definition: SpecialPixel.h:199
Obtain SPICE rotation information for a body.
This class is used to modify and manage solve settings for 1 to many BundleObservations.
QString m_instrumentId
Spacecraft instrument id.
Isis exception class.
Definition: IException.h:99
InstrumentPositionSolveOption
Options for how to solve for instrument position.
Obtain SPICE position information for a body.
SpicePosition * m_instrumentPosition
Instrument spice position (in primary image).
InstrumentPointingSolveOption
Options for how to solve for instrument pointing.
double GetTimeScale()
Accessor method to get the rotation time scale.
void SetPolynomial(const Source type=PolyFunction)
Set the coefficients of a polynomial fit to each of the three camera angles for the time period cover...
const double RAD2DEG(57.29577951308232087679815481)
Multiplier for converting from radians to degrees.
int m_index
Index of this observation.
void SetPolynomial(const Source type=PolyFunction)
Set the coefficients of a polynomial fit to each of the components (X, Y, Z) of the position vector f...

U.S. Department of the Interior | U.S. Geological Survey
ISIS | Privacy & Disclaimers | Astrogeology Research Program
To contact us, please post comments and questions on the ISIS Support Center
File Modified: 07/12/2023 23:15:06