GeneralBrokenLines V03-01-03
using EIGEN
GblData.cpp
Go to the documentation of this file.
1/*
2 * GblData.cpp
3 *
4 * Created on: Aug 18, 2011
5 * Author: kleinwrt
6 */
7
30#include "GblData.h"
31using namespace Eigen;
32
34namespace gbl {
35
37
46GblData::GblData(unsigned int aLabel, dataBlockType aType, double aValue,
47 double aPrec, unsigned int aTraj, unsigned int aPoint,
48 unsigned int aMeas) :
49 theLabel(aLabel), theRow(0), theType(aType), theValue(aValue), thePrecision(
50 aPrec), theTrajectory(aTraj), thePoint(aPoint), theMeas(aMeas), theDWMethod(
51 0), theDownWeight(1.), thePrediction(0.), theNumLocal(0), moreParameters(), moreDerivatives() {
52
53}
54
56}
57
59
64void GblData::addDerivatives(const std::vector<unsigned int> &index,
65 const std::vector<double> &derivatives) {
66 for (unsigned int i = 0; i < derivatives.size(); ++i) // any derivatives
67 {
68 if (derivatives[i]) {
69 moreParameters.push_back(index[i]);
70 moreDerivatives.push_back(derivatives[i]);
71 }
72 }
73}
74
76void GblData::setPrediction(const VVector &aVector) {
77
78 thePrediction = 0.;
79 if (theNumLocal > 0) {
80 for (unsigned int i = 0; i < theNumLocal; ++i) {
81 thePrediction += theDerivatives[i] * aVector(theParameters[i] - 1);
82 }
83 } else {
84 for (unsigned int i = 0; i < moreDerivatives.size(); ++i) {
86 * aVector(moreParameters[i] - 1);
87 }
88 }
89}
90
92
95double GblData::setDownWeighting(unsigned int aMethod) {
96
97 theDWMethod = aMethod;
98 double aWeight = 1.;
99 double scaledResidual = fabs(theValue - thePrediction) * sqrt(thePrecision);
100 if (aMethod == 1) // Tukey
101 {
102 if (scaledResidual < 4.6851) {
103 aWeight = (1.0 - 0.045558 * scaledResidual * scaledResidual);
104 aWeight *= aWeight;
105 } else {
106 aWeight = 0.;
107 }
108 } else if (aMethod == 2) //Huber
109 {
110 if (scaledResidual >= 1.345) {
111 aWeight = 1.345 / scaledResidual;
112 }
113 } else if (aMethod == 3) //Cauchy
114 {
115 aWeight = 1.0 / (1.0 + (scaledResidual * scaledResidual / 5.6877));
116 }
117 theDownWeight = aWeight;
118 return aWeight;
119}
120
122
127double GblData::getChi2() const {
128 double scaledResidual = fabs(theValue - thePrediction) * sqrt(thePrecision);
129 double chi2 = scaledResidual * scaledResidual;
130 if (theDWMethod == 1) // Tukey
131 {
132 if (scaledResidual < 4.6851) {
133 chi2 = (1.0
134 - pow(1.0 - 0.045558 * scaledResidual * scaledResidual, 3))
135 / (3. * 0.045558);
136 } else {
137 chi2 = 1.0 / (3. * 0.045558);
138 }
139 } else if (theDWMethod == 2) //Huber
140 {
141 if (scaledResidual >= 1.345) {
142 chi2 = 1.345 * (2. * scaledResidual - 1.345);
143 }
144 } else if (theDWMethod == 3) //Cauchy
145 {
146 chi2 = log(1.0 + (scaledResidual * scaledResidual / 5.6877)) * 5.6877;
147 }
148 return chi2;
149}
150
152void GblData::printData() const {
153
155 std::cout << " measurement at label " << theLabel << " of type "
156 << theType << " from meas, row " << theMeas << ", " << theRow
157 << ": " << theValue << ", " << thePrecision << std::endl;
158 } else {
159 std::cout << " measurement at label " << theLabel << " of type "
160 << theType << " from row " << theRow << ": " << theValue << ", "
161 << thePrecision << std::endl;
162 }
163 std::cout << " param " << moreParameters.size() + theNumLocal << ":";
164 for (unsigned int i = 0; i < moreParameters.size(); ++i) {
165 std::cout << " " << moreParameters[i];
166 }
167 for (unsigned int i = 0; i < theNumLocal; ++i) {
168 std::cout << " " << theParameters[i];
169 }
170 std::cout << std::endl;
171 std::cout << " deriv " << moreDerivatives.size() + theNumLocal << ":";
172 for (unsigned int i = 0; i < moreDerivatives.size(); ++i) {
173 std::cout << " " << moreDerivatives[i];
174 }
175 for (unsigned int i = 0; i < theNumLocal; ++i) {
176 std::cout << " " << theDerivatives[i];
177 }
178 std::cout << std::endl;
179}
180
182
185unsigned int GblData::getLabel() const {
186 return theLabel;
187}
188
190
194 return theType;
195}
196
198
205void GblData::getLocalData(double &aValue, double &aWeight,
206 unsigned int &numLocal, unsigned int *&indLocal, double *&derLocal) {
207
208 aValue = theValue;
209 aWeight = thePrecision * theDownWeight;
210 if (theNumLocal > 0) {
211 numLocal = theNumLocal;
212 indLocal = theParameters;
213 derLocal = theDerivatives;
214 } else {
215 numLocal = moreParameters.size();
216 indLocal = &moreParameters[0];
217 derLocal = &moreDerivatives[0];
218 }
219}
220
222
233void GblData::getAllData(double &aValue, double &aErr, unsigned int &numLocal,
234 unsigned int *&indLocal, double *&derLocal, unsigned int &aTraj,
235 unsigned int &aPoint, unsigned int &aMeas, unsigned int &aRow) {
236 aValue = theValue;
237 aErr = 1.0 / sqrt(thePrecision);
238 if (theNumLocal > 0) {
239 numLocal = theNumLocal;
240 indLocal = theParameters;
241 derLocal = theDerivatives;
242 } else {
243 numLocal = moreParameters.size();
244 indLocal = &moreParameters[0];
245 derLocal = &moreDerivatives[0];
246 }
247 aTraj = theTrajectory;
248 aPoint = thePoint;
249 aMeas = theMeas;
250 aRow = theRow;
251}
252
254
262void GblData::getResidual(double &aResidual, double &aVariance,
263 double &aDownWeight, unsigned int &numLocal, unsigned int *&indLocal,
264 double *&derLocal) {
265 aResidual = theValue - thePrediction;
266 aVariance = 1.0 / thePrecision;
267 aDownWeight = theDownWeight;
268 if (theNumLocal > 0) {
269 numLocal = theNumLocal;
270 indLocal = theParameters;
271 derLocal = theDerivatives;
272 } else {
273 numLocal = moreParameters.size();
274 indLocal = &moreParameters[0];
275 derLocal = &moreDerivatives[0];
276 }
277}
278
280
284void GblData::getResidual(double &aResidual, double &aVariance) {
285 aResidual = theValue - thePrediction;
286 aVariance = 1.0 / thePrecision;
287}
288}
GblData definition.
double thePrecision
Precision (1/sigma**2)
Definition: GblData.h:130
std::vector< unsigned int > moreParameters
List of fit parameters (with non zero derivatives)
Definition: GblData.h:142
unsigned int theMeas
Measurement number (at point)
Definition: GblData.h:133
double theDerivatives[9]
List of derivatives for fit.
Definition: GblData.h:140
double theDownWeight
Down-weighting factor (0-1)
Definition: GblData.h:135
double thePrediction
Prediction from fit.
Definition: GblData.h:136
unsigned int theTrajectory
Trajectory number.
Definition: GblData.h:131
unsigned int theNumLocal
Number of (non zero) local derivatives (max 9 for kinks+steps)
Definition: GblData.h:138
unsigned int theLabel
Label (of corresponding point)
Definition: GblData.h:126
dataBlockType theType
Type (None, InternalMeasurement, InternalKink, ExternalSeed, ExternalMeasurement)
Definition: GblData.h:128
unsigned int theDWMethod
Down-weighting method (0: None, 1: Tukey, 2: Huber, 3: Cauchy)
Definition: GblData.h:134
virtual ~GblData()
Definition: GblData.cpp:55
void addDerivatives(unsigned int iRow, const std::array< unsigned int, 5 > &labDer, const Matrix5d &matDer, unsigned int iOff, const Eigen::MatrixBase< LocalDerivative > &derLocal, unsigned int extOff, const Eigen::MatrixBase< ExtDerivative > &extDer)
Add derivatives from measurement.
Definition: GblData.h:147
dataBlockType getType() const
Get type.
Definition: GblData.cpp:193
unsigned int getLabel() const
Get label.
Definition: GblData.cpp:185
unsigned int theRow
Row number (of measurement)
Definition: GblData.h:127
double getChi2() const
Calculate Chi2 contribution.
Definition: GblData.cpp:127
void printData() const
Print data block.
Definition: GblData.cpp:152
void getResidual(double &aResidual, double &aVariance, double &aDownWeight, unsigned int &numLocal, unsigned int *&indLocal, double *&derLocal)
Get data for residual (and errors) "long list".
Definition: GblData.cpp:262
unsigned int thePoint
Point number (on trajectory)
Definition: GblData.h:132
GblData(unsigned int aLabel, dataBlockType aType, double aValue, double aPrec, unsigned int aTraj=0, unsigned int aPoint=0, unsigned int aMeas=0)
Create data block.
Definition: GblData.cpp:46
void getLocalData(double &aValue, double &aWeight, unsigned int &numLocal, unsigned int *&indLocal, double *&derLocal)
Get Data for local fit.
Definition: GblData.cpp:205
void getAllData(double &aValue, double &aErr, unsigned int &numLocal, unsigned int *&indLocal, double *&derLocal, unsigned int &aTraj, unsigned int &aPoint, unsigned int &aMeas, unsigned int &aRow)
Get all Data for MP-II binary record.
Definition: GblData.cpp:233
double theValue
Value (residual)
Definition: GblData.h:129
unsigned int theParameters[9]
List of parameters (with non zero derivatives)
Definition: GblData.h:139
double setDownWeighting(unsigned int aMethod)
Outlier down weighting with M-estimators (by GblTrajectory::fit).
Definition: GblData.cpp:95
std::vector< double > moreDerivatives
List of derivatives for fit.
Definition: GblData.h:143
void setPrediction(const VVector &aVector)
Calculate prediction for data from fit (by GblTrajectory::fit).
Definition: GblData.cpp:76
Simple Vector based on std::vector<double>
Definition: VMatrix.h:43
Namespace for the general broken lines package.
dataBlockType
Definition: GblData.h:49
@ InternalMeasurement
Definition: GblData.h:50