Isis Developer Reference
ZeroBufferFit.h
Go to the documentation of this file.
1 #ifndef ZeroBufferFit_h
2 #define ZeroBufferFit_h
3 
10 /* SPDX-License-Identifier: CC0-1.0 */
11 
12 #include <string>
13 #include <vector>
14 
15 #include "HiCalTypes.h"
16 #include "HiCalUtil.h"
17 #include "HiCalConf.h"
18 #include "NonLinearLSQ.h"
19 #include "Module.h"
20 #include "IException.h"
21 
22 namespace Isis {
23 
42  class ZeroBufferFit : public NonLinearLSQ, public Module {
43 
44  public:
45  ZeroBufferFit(const HiCalConf &conf);
47  virtual ~ZeroBufferFit() { }
48 
54  void setBin(int bin) { _timet.setBin(bin); }
55 
61  void setLineTime(double ltime) { _timet.setLineTime(ltime); }
62 
71  inline int size() const { return (_data.dim()); }
72 
84  int nSize() const { return (_b2.dim()); }
85 
94  int nParms() const { return (4); }
95 
97  void setabsErr(double absError) { _absErr = absError; }
99  void setrelErr(double relError) { _relErr = relError; }
101  double absErr() const { return (_absErr); }
103  double relErr() const { return (_relErr); }
104 
105  HiVector Solve(const HiVector &d);
106  NLVector guess();
107  int checkIteration(const int Iter, const NLVector &fitcoefs,
108  const NLVector &uncerts, double cplxconj,
109  int Istatus);
110 
111  NLVector f_x(const NLVector &a);
112  NLMatrix df_x(const NLVector &a);
113 
115  double Chisq() const { return (_chisq); }
117  int DoF() const { return (nSize() - nParms()); }
118 
119  HiVector Yfit() const;
120  HiVector Normalize(const HiVector &v);
121 
122  private:
123  HiLineTimeEqn _timet; // This is the X data set
124  HiVector _data; // Typically will be the HiRISE buffer data
125  HiVector _b2; // Data buffer used in fitting
126  double _absErr; // Absolute error convergence test
127  double _relErr; // Relative error convergence test
128  double _maxLog; // Maximum log value to constrain
129  int _badLines; // Exclude lines at end of buffer
130  int _sWidth; // Width of guestimate filter
131  int _sIters; // Filter interations
132  bool _skipFit; // Skip fitting and pass input through
133  bool _useLinFit; // Use linear fit on failure of LM, else Zf
134  int _minLines; // Minimum number of lines to fit (Default: 100)
135  HiVector _cc; // Parameter of 2-D fit
136  HiVector _guess; // Initial guestimate of solutions
137  HiVector _coefs; // Coefficients of solution
138  HiVector _uncert; // Uncertanties
139  double _chisq; // ChiSq of NonLinear equation
140 
141  HiVector poly_fit(const HiVector &d, const double line0 = 0.0) const;
142  virtual void printOn(std::ostream &o) const;
144  int goodLines(const HiVector &d) const { return (d.dim() - _badLines); }
146  bool gotGoodLines(const HiVector &d) const {
147  return (goodLines(d) >= _minLines);
148  }
149  };
150 
151 } // namespace Isis
152 #endif
Isis::Module::_history
HiHistory _history
Hierarchial component history.
Definition: Module.h:152
Isis::ZeroBufferFit::setrelErr
void setrelErr(double relError)
Sets the relative error parameter
Definition: ZeroBufferFit.h:99
Isis::IsEqual
bool IsEqual(const QString &v1, const QString &v2="TRUE")
Shortened string equality test.
Definition: HiCalUtil.h:258
Isis::Module::_fmtWidth
int _fmtWidth
Default field with of double.
Definition: Module.h:153
Isis::Statistics
This class is used to accumulate statistics on double arrays.
Definition: Statistics.h:94
Isis::Statistics::AddData
void AddData(const double *data, const unsigned int count)
Add an array of doubles to the accumulators and counters.
Definition: Statistics.cpp:141
Isis::DbProfile::Name
QString Name() const
Returns the name of this property.
Definition: DbProfile.h:104
Isis::ZeroBufferFit::Solve
HiVector Solve(const HiVector &d)
Compute non-linear fit to (typically) ZeroBufferSmooth module.
Definition: ZeroBufferFit.cpp:91
Isis::ToString
QString ToString(const T &value)
Helper function to convert values to strings.
Definition: HiCalUtil.h:246
Isis::DbProfile::exists
bool exists(const QString &key) const
Checks for the existance of a keyword.
Definition: DbProfile.h:115
Isis::ZeroBufferFit::size
int size() const
Returns the size of the data buffer.
Definition: ZeroBufferFit.h:71
Isis::Module
Module manages HiRISE calibration vectors from various sources.
Definition: Module.h:39
Isis::Module::ref
const HiVector & ref() const
Return data via a const reference.
Definition: Module.h:111
Isis::ConfKey
T ConfKey(const DbProfile &conf, const QString &keyname, const T &defval, int index=0)
Find a keyword in a profile using default for non-existant keywords.
Definition: HiCalUtil.h:205
Isis::ZeroBufferFit::ZeroBufferFit
ZeroBufferFit(const HiCalConf &conf)
Compute second level drift correction (Zf module)
Definition: ZeroBufferFit.cpp:38
Isis::ZeroBufferFit::setabsErr
void setabsErr(double absError)
Sets the absolute error parameter.
Definition: ZeroBufferFit.h:97
HiCalConf.h
Isis::ZeroBufferFit
Computes non-linear lsq fit of HiRISE Drift (Zd module)
Definition: ZeroBufferFit.h:42
Isis::Statistics::Reset
void Reset()
Reset all accumulators and counters to zero.
Definition: Statistics.cpp:113
Isis::NonLinearLSQ::curvefit
int curvefit()
Definition: NonLinearLSQ.cpp:28
Isis::ZeroBufferFit::absErr
double absErr() const
Returns the current value of the absolute error.
Definition: ZeroBufferFit.h:101
Isis::NonLinearLSQ::uncert
NLVector uncert() const
Return uncertainties from last fit processing.
Definition: NonLinearLSQ.h:105
Isis::ZeroBufferFit::setLineTime
void setLineTime(double ltime)
Set scan line time.
Definition: ZeroBufferFit.h:61
Isis::NonLinearLSQ::setMaxIters
void setMaxIters(int m)
Sets the maximum number of iterations.
Definition: NonLinearLSQ.h:62
ZeroBufferFit.h
IString.h
Isis::LowPassFilter
Compute a low pass filter from a Module class content.
Definition: LowPassFilter.h:30
Isis::MultivariateStatistics::AddData
void AddData(const double *x, const double *y, const unsigned int count)
Add two arrays of doubles to the accumulators and counters.
Definition: MultivariateStatistics.cpp:66
Isis::NonLinearLSQ::coefs
NLVector coefs() const
Return coefficients from last fit processing.
Definition: NonLinearLSQ.h:103
LowPassFilter.h
Isis::ZeroBufferFit::DoF
int DoF() const
Returns the Degrees of Freedom.
Definition: ZeroBufferFit.h:117
Isis::NonLinearLSQ::NLVector
TNT::Array1D< double > NLVector
Definition: NonLinearLSQ.h:44
Isis::Module::formatDbl
QString formatDbl(const double &value) const
Properly format values that could be special pixels.
Definition: Module.h:169
Isis::ZeroBufferFit::Normalize
HiVector Normalize(const HiVector &v)
Compute normalized solution vector from result.
Definition: ZeroBufferFit.cpp:282
Isis::NonLinearLSQ::success
bool success() const
Determine success from last fit processing.
Definition: NonLinearLSQ.h:83
Isis::NonLinearLSQ
NonLinearLSQ Computes a fit using a Levenberg-Marquardt algorithm.
Definition: NonLinearLSQ.h:42
Isis::MultivariateStatistics
Container of multivariate statistics.
Definition: MultivariateStatistics.h:54
Isis::toInt
int toInt(const QString &string)
Global function to convert from a string to an integer.
Definition: IString.cpp:93
Isis::NonLinearLSQ::maxIters
int maxIters() const
Maximum number iterations for valid solution.
Definition: NonLinearLSQ.h:70
Isis::HiLineTimeEqn::setLineTime
void setLineTime(double ltime)
Definition: HiCalUtil.h:367
Isis::HiLineTimeEqn
Compute HiRISE line times.
Definition: HiCalUtil.h:361
Isis::ToInteger
int ToInteger(const T &value)
Helper function to convert values to Integers.
Definition: HiCalUtil.h:222
Isis::ZeroBufferFit::df_x
NLMatrix df_x(const NLVector &a)
Computes the first derivative of the function at the current iteration.
Definition: ZeroBufferFit.cpp:245
Isis::IsTrueValue
bool IsTrueValue(const DbProfile &prof, const QString &key, const QString &value="TRUE")
Determines if the keyword value is the expected value.
Definition: HiCalUtil.h:276
HiCalUtil.h
Isis::MIN
T MIN(const T &A, const T &B)
Implement templatized MIN fumnction.
Definition: PhotometricFunction.h:34
LeastSquares.h
Isis::DbProfile
A DbProfile is a container for access parameters to a database.
Definition: DbProfile.h:51
HiCalTypes.h
Isis::HiVector
TNT::Array1D< double > HiVector
1-D Buffer
Definition: HiCalTypes.h:27
NonLinearLSQ.h
Isis::HiHistory::add
void add(const QString &event)
Definition: HiCalTypes.h:55
Isis::NonLinearLSQ::statusstr
std::string statusstr() const
Return error message pertaining to last fit procesing.
Definition: NonLinearLSQ.h:87
MultivariateStatistics.h
Isis::ZeroBufferFit::setBin
void setBin(int bin)
Set binning/summing mode.
Definition: ZeroBufferFit.h:54
Isis::ZeroBufferFit::nParms
int nParms() const
Number of parameter to be fitted.
Definition: ZeroBufferFit.h:94
Isis::Statistics::Average
double Average() const
Computes and returns the average.
Definition: Statistics.cpp:300
Isis::ToDouble
double ToDouble(const T &value)
Helper function to convert values to doubles.
Definition: HiCalUtil.h:234
Isis::MultivariateStatistics::LinearRegression
void LinearRegression(double &a, double &b) const
Fits a line.
Definition: MultivariateStatistics.cpp:222
Isis::ZeroBufferFit::Yfit
HiVector Yfit() const
Computes the solution vector using current coefficents.
Definition: ZeroBufferFit.cpp:266
Isis::ZeroBufferFit::relErr
double relErr() const
Returns the current value of the relative error.
Definition: ZeroBufferFit.h:103
IException.h
Isis::toDouble
double toDouble(const QString &string)
Global function to convert from a string to a double.
Definition: IString.cpp:149
std
Namespace for the standard library.
Isis::ZeroBufferFit::checkIteration
int checkIteration(const int Iter, const NLVector &fitcoefs, const NLVector &uncerts, double cplxconj, int Istatus)
Computes the interation check for convergence.
Definition: ZeroBufferFit.cpp:219
Isis::NonLinearLSQ::NLMatrix
TNT::Array2D< double > NLMatrix
Definition: NonLinearLSQ.h:45
Isis::ZeroBufferFit::nSize
int nSize() const
Returns the size of the fitted buffer.
Definition: ZeroBufferFit.h:84
Isis::HiLineTimeEqn::setBin
void setBin(int bin)
Definition: HiCalUtil.h:368
Isis::ZeroBufferFit::f_x
NLVector f_x(const NLVector &a)
Computes the function value at the current iteration.
Definition: ZeroBufferFit.cpp:227
Isis::ZeroBufferFit::guess
NLVector guess()
Compute the initial guess of the fit.
Definition: ZeroBufferFit.cpp:162
Isis::ZeroBufferFit::~ZeroBufferFit
virtual ~ZeroBufferFit()
Destructor.
Definition: ZeroBufferFit.h:47
Module.h
Isis::NonLinearLSQ::nIterations
int nIterations() const
Return number of iterations from last fit processing.
Definition: NonLinearLSQ.h:107
Isis::ZeroBufferFit::Chisq
double Chisq() const
Returns the Chi-Square value of the fit solution.
Definition: ZeroBufferFit.h:115
Isis
This is free and unencumbered software released into the public domain.
Definition: Apollo.h:16