Commit 6ffeb278 authored by dauncey's avatar dauncey

Missing files

parent 3269c598
Pipeline #14676 skipped
No preview for this file type
......@@ -61,6 +61,23 @@ Inverting the large matrix $M$ is required to give the solution for the
optimal $a_i$.
Note, $M$ is similar (but not identical) to the error matrix of the $d_i$.
The resulting RMS using the best fit values is given by
\begin{eqnarray*}
\mathrm{RMS}^2_\mathrm{min}
&=& \frac{1}{N} \sum_e \left[
\left(\sum_i a_i d_{e,i} \right)^2
- 2 T_e \sum_i a_i d_{e,i} + T_e^2 \right] \\
&=& \frac{1}{N} \sum_j \sum_i a_j a_i \sum_e d_{e,j} d_{e,i}
- \frac{2}{N} \sum_i a_i \sum_e T_e d_{e,i} + \frac{1}{N} \sum_e T_e^2 \\
&=& \sum_j \sum_i a_j a_i M_{ji}
- 2 \sum_i a_i v_i + \frac{1}{N} \sum_e T_e^2
= a^T M a - 2 a^T v + \frac{1}{N} \sum_e T_e^2
\end{eqnarray*}
But since the solution is defined by $Ma=v$, then $a^T M a = a^T v$. Hence
\begin{equation}
\mathrm{RMS}^2_\mathrm{min} = \frac{1}{N} \sum_e T_e^2 - a^T M a
= \frac{1}{N} \sum_e T_e^2 - v^T M^{-1} v
\end{equation}
\section{Units}
Keeping quantities to 16-bit integers.
......
......@@ -323,6 +323,50 @@ public:
200,0.0,200.0);
}
// SimHits; all 1k events fitted with 200 PU
fFitCoeff[0]=-2.24876;
fFitCoeff[1]=1.28108;
fFitCoeff[2]=1.01436;
fFitCoeff[3]=2.49601;
fFitCoeff[4]=1.98179;
fFitCoeff[5]=1.20344;
fFitCoeff[6]=1.00278;
fFitCoeff[7]=0.695715;
fFitCoeff[8]=1.19054;
fFitCoeff[9]=1.68355;
fFitCoeff[10]=0.800412;
fFitCoeff[11]=1.32589;
fFitCoeff[12]=0.91138;
fFitCoeff[13]=1.3035;
fFitCoeff[14]=0.689325;
fFitCoeff[15]=1.25186;
fFitCoeff[16]=1.46017;
fFitCoeff[17]=1.68828;
fFitCoeff[18]=0.939856;
fFitCoeff[19]=0.682892;
fFitCoeff[20]=0.434451;
fFitCoeff[21]=1.63003;
fFitCoeff[22]=1.94675;
fFitCoeff[23]=2.39824;
fFitCoeff[24]=1.11453;
fFitCoeff[25]=-0.842631;
fFitCoeff[26]=-0.687827;
fFitCoeff[27]=-0.360403;
fFitCoeff[28]=0;
fFitCoeff[29]=0;
fFitCoeff[30]=0;
fFitCoeff[31]=0;
fFitCoeff[32]=0;
fFitCoeff[33]=0;
fFitCoeff[34]=0;
fFitCoeff[35]=0;
fFitCoeff[36]=0;
fFitCoeff[37]=0;
fFitCoeff[38]=0;
fFitCoeff[39]=0;
// SimHits: all 1k events in SinglePhotonPt35_simHit_1kev
fSgShape[0]=0.00261114;
fSgShape[1]=0.00536047;
......@@ -2393,7 +2437,7 @@ public:
bool fitShape(const double *a, double &et) {
unsigned method(1);
unsigned method(2);
if(method==0) {
double chiMin(1.0e100);
......@@ -2495,6 +2539,15 @@ public:
<< " and PU = " << pu << std::endl;
}
if(method==2) {
et=0.0;
for(unsigned i(0);i<40;i++) {
et+=a[i]*fFitCoeff[i];
}
std::cout << "Fit coefficients give transverse energy = " << et
<< std::endl;
}
return true;
}
......@@ -2599,8 +2652,8 @@ public:
double tanMax[2];
double phiMax[2];
simEvent(event);
return true;
//simEvent(event);
//return true;
for(unsigned e(0);e<Geometry::kNumberOfEndcaps;e++) {
......@@ -3157,6 +3210,7 @@ protected:
TProfile *pDeltaEVsRho[40],*pDeltaE2VsRho[40],
*pSelDeltaEVsRho[40],*pSelDeltaE2VsRho[40];
double fFitCoeff[40];
double fSgShape[40];
double fSdShape[40];
double fPuShape[40];
......
#ifndef ErrorMatrix_HH
#define ErrorMatrix_HH
#include <iostream>
#include <iomanip>
#include <fstream>
#include <cstdint>
#include <cassert>
#include <cmath>
#include "TMatrixDSym.h"
class ErrorMatrix {
public:
ErrorMatrix() {
fNum=0.0;
for(unsigned i(0);i<40;i++) {
fSum[i]=0.0;
for(unsigned j(0);j<40;j++) {
fSSq[i][j]=0.0;
}
}
}
void fill(const double *p) {
fNum++;
for(unsigned i(0);i<40;i++) {
fSum[i]+=p[i];
for(unsigned j(0);j<40;j++) {
fSSq[i][j]+=p[i]*p[j];
}
}
}
void normalise() {
assert(fNum>0.5);
for(unsigned i(0);i<40;i++) {
fSum[i]/=fNum;
}
for(unsigned i(0);i<40;i++) {
for(unsigned j(0);j<40;j++) {
fSSq[i][j]/=fNum;
}
}
fNum=1.0;
}
TMatrixDSym fitMatrix() const {
TMatrixDSym fm(40);
for(unsigned i(0);i<40;i++) {
for(unsigned j(0);j<40;j++) {
fm(i,j)=fSSq[i][j];
}
}
return fm;
}
TMatrixDSym errorMatrix() const {
TMatrixDSym em(40);
for(unsigned i(0);i<40;i++) {
for(unsigned j(0);j<40;j++) {
em(i,j)=fSSq[i][j]-fSum[i]*fSum[j];
}
}
return em;
}
TMatrixDSym weightMatrix() const {
TMatrixDSym wm(errorMatrix());
wm.Invert();
return wm;
}
void read(const std::string &fName) {
std::ifstream fin(fName.c_str());
if(!fin) return;
fin >> fNum;
for(unsigned i(0);i<40;i++) {
fin >> fSum[i];
}
for(unsigned i(0);i<40;i++) {
for(unsigned j(0);j<40;j++) {
fin >> fSSq[i][j];
}
}
}
void write(const std::string &fName) {
std::ofstream fout(fName.c_str());
if(!fout) return;
fout << fNum << std::endl;
for(unsigned i(0);i<40;i++) {
fout << " " << std::setprecision(12) << fSum[i];
}
fout << std::endl;
for(unsigned i(0);i<40;i++) {
for(unsigned j(0);j<40;j++) {
fout << " " << std::setprecision(12) << fSSq[i][j];
}
fout << std::endl;
}
}
double fNum;
double fSum[40];
double fSSq[40][40];
};
#endif
#include <iostream>
#include <sstream>
#include "TObject.h"
#include "ErrorMatrix.hh"
//TMatrix x;
int main(int argc, char* argv[]) {
//TMatrixDSym sm;
ErrorMatrix em;
em.read("AnalysisPhoton_ErrorMatrix_Save.out");
em.normalise();
//em.weightMatrix();
return 0;
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment