diff --git a/NPLib/Calibration/NPSiliconCalibrator.cxx b/NPLib/Calibration/NPSiliconCalibrator.cxx
index 1242e5b63659b2eb7618ceaa60b11b631d3afb43..36459217c9e9d99b14757214e8fb0d98c96d697a 100644
--- a/NPLib/Calibration/NPSiliconCalibrator.cxx
+++ b/NPLib/Calibration/NPSiliconCalibrator.cxx
@@ -24,7 +24,7 @@ NPL::SiliconCalibrator::~SiliconCalibrator(){
}
//////////////////////////////////////////////
-double NPL::SiliconCalibrator::ZeroExtrapolation(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, unsigned int pedestal, unsigned int max_iteration,double rmin,double rmax){
+double NPL::SiliconCalibrator::ZeroExtrapolation(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, unsigned int pedestal, unsigned int max_iteration, double rmin,double rmax){
if(histo->GetEntries()==0){
coeff.clear();
coeff.push_back(0);
@@ -80,8 +80,12 @@ double NPL::SiliconCalibrator::ZeroExtrapolation(TH1* histo, NPL::CalibrationSou
}
Assume_E.resize(source_size,0);
TGraphErrors* g = FitSpectrum(histo,rmin,rmax);
- if (g->GetN()!=3)
+ if (g->GetN()!=3){
+ coeff.clear();
+ coeff.push_back(0);
+ coeff.push_back(-1);
return -1;
+ }
while(k++<max_iteration && abs(Assume_Thickness) < 100*micrometer){
@@ -92,7 +96,7 @@ double NPL::SiliconCalibrator::ZeroExtrapolation(TH1* histo, NPL::CalibrationSou
g->SetPoint(i,x[i] ,Assume_E[i]);
}
- DistanceToPedestal = FitPoints(g,&Assume_E[0] , Source_Sig, coeff, 0 );
+ DistanceToPedestal = FitPoints(g,&Assume_E[0], Source_Sig, coeff, 0 );
if(abs(DistanceToPedestal) < my_precision || Step < 0.01*micrometer)
break;
@@ -134,7 +138,7 @@ double NPL::SiliconCalibrator::ZeroExtrapolation(TH1* histo, NPL::CalibrationSou
return Assume_Thickness/micrometer;
}
////////////////////////////////////////////////////////////////////////////////
-double NPL::SiliconCalibrator::SimpleCalibration(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, double rmin,double rmax){
+double NPL::SiliconCalibrator::SimpleCalibration(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, double Assume_Thickness, double rmin,double rmax){
if(histo->GetEntries()==0){
coeff.clear();
coeff.push_back(0);
@@ -157,13 +161,14 @@ double NPL::SiliconCalibrator::SimpleCalibration(TH1* histo, NPL::CalibrationSou
}
}
- if(counts < 30){
+ if(counts == 0){
coeff.clear();
coeff.push_back(0);
coeff.push_back(-1);
return -2;
}
+
m_CalibrationSource= CS;
m_EL_Al= EL;
@@ -175,14 +180,31 @@ double NPL::SiliconCalibrator::SimpleCalibration(TH1* histo, NPL::CalibrationSou
Source_E[i] = CS->GetEnergies()[i][0];
Source_Sig[i] = CS->GetEnergiesErrors()[i][0];
}
- TGraphErrors* g = FitSpectrum(histo,rmin,rmax);
- double dist = FitPoints(g,Source_E , Source_Sig, coeff, 0 );
+
+ TGraphErrors* g = FitSpectrum(histo,rmin,rmax); // the graph uses the original tabulated alpha energies
+ if(g->GetN() != 3) {
+ coeff.clear();
+ coeff.push_back(0);
+ coeff.push_back(-1);
+ return -1;
+ }
+
+ // Compute the new assumed energies and store in graph
+ double* x = g->GetX(); // get the fitted channel centroids
+ vector<double> Assume_E ; // Energie calculated assuming Assume_Thickness deadlayer of Al
+ Assume_E.resize(source_size,0);
+ for(unsigned int i = 0 ; i < source_size ; i++){
+ Assume_E[i] = EL->Slow(Source_E[i], Assume_Thickness, 0); // Assume_Thickness in micrometer
+ g->SetPoint(i,x[i],Assume_E[i]);
+ }
+
+ double DistanceToPedestal = FitPoints(g, &Assume_E[0] , Source_Sig, coeff, 0 );
delete g;
- if(dist!=-3)
- return abs(dist);
+ if(DistanceToPedestal!=-3)
+ return abs(DistanceToPedestal);
else
- return dist;
+ return DistanceToPedestal;
}
@@ -196,7 +218,7 @@ double NPL::SiliconCalibrator::FitPoints(TGraphErrors* gre, double* Energies, do
coeff.push_back(gre->GetFunction("pol1")->GetParameter(0));
coeff.push_back(gre->GetFunction("pol1")->GetParameter(1));
// Compute the Distance to pedestal:
- return (coeff[0]/coeff[1]-pedestal);
+ return ((-coeff[0]/coeff[1])-pedestal);
}
else {
//std::cout << "gre->GetN() <= 0" << std::endl; // used for debugging
@@ -281,8 +303,9 @@ TGraphErrors* NPL::SiliconCalibrator::FitSpectrum(TH1* histo, double rmin, doubl
return gre;
}
- else
+ else{
return (new TGraphErrors());
+ }
}
diff --git a/NPLib/Calibration/NPSiliconCalibrator.h b/NPLib/Calibration/NPSiliconCalibrator.h
index cdc5b97a99852be0362283c72ad2bffe9edb5ce2..eab6a935d49b4d1cd46bbb3efab47da19fdfeb03 100644
--- a/NPLib/Calibration/NPSiliconCalibrator.h
+++ b/NPLib/Calibration/NPSiliconCalibrator.h
@@ -26,7 +26,7 @@ namespace NPL{
// Use the Zero Extrapolation method to perform fit and return the dead layer thickness
double ZeroExtrapolation(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, unsigned int pedestal, unsigned int max_iteration = 10000 , double rmin=-1,double rmax=-1);
- double SimpleCalibration(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, double rmin=-1,double rmax=-1);
+ double SimpleCalibration(TH1* histo, NPL::CalibrationSource* CS, NPL::EnergyLoss* EL, vector<double>& coeff, double Assume_Thickness=0,double rmin=-1,double rmax=-1);
// Return distance to pedestal. Use energies in Energies to perform fit and fill coeff with the results
double FitPoints(TGraphErrors* Graph, double* Energies , double* ErrEnergies, vector<double>& coeff , double pedestal = 0 );
diff --git a/NPLib/Detectors/GeTAMU/TGeTAMUPhysics.cxx b/NPLib/Detectors/GeTAMU/TGeTAMUPhysics.cxx
index 70350ca756dcef4829fe3b52f249a55b15f7f597..e56285acb5d2896897b90be8706167ba50a19aa1 100644
--- a/NPLib/Detectors/GeTAMU/TGeTAMUPhysics.cxx
+++ b/NPLib/Detectors/GeTAMU/TGeTAMUPhysics.cxx
@@ -38,6 +38,9 @@ using namespace NPUNITS;
// ROOT
#include "TChain.h"
+#include "TRandom3.h"
+
+TRandom *Random = new TRandom3();
///////////////////////////////////////////////////////////////////////////
ClassImp(TGeTAMUPhysics)
@@ -194,7 +197,7 @@ void TGeTAMUPhysics::PreTreat(){
clover = m_EventData->GetCoreCloverNbrE(i);
crystal = m_EventData->GetCoreCrystalNbrE(i);
name = "GETAMU/D"+ NPL::itoa(clover)+"_CRY"+ NPL::itoa(crystal);
- Energy = cal->ApplyCalibration(name+"_E", Eraw);
+ Energy = cal->ApplyCalibration(name+"_E", Eraw+Random->Rndm());
Singles_CloverMap_CryEN[clover].push_back(crystal);
Singles_CloverMap_CryE[clover].push_back(Energy);
m_PreTreatedData->SetCoreE(clover,crystal,Energy);
@@ -207,7 +210,7 @@ void TGeTAMUPhysics::PreTreat(){
clover = m_EventData->GetCoreCloverNbrT(i);
crystal = m_EventData->GetCoreCrystalNbrT(i);
name = "GETAMU/D"+ NPL::itoa(clover)+"_CRY"+ NPL::itoa(crystal);
- Time = cal->ApplyCalibration(name+"_T", Traw);
+ Time = cal->ApplyCalibration(name+"_T", Traw+Random->Rndm());
Singles_CloverMap_CryTN[clover].push_back(crystal);
Singles_CloverMap_CryT[clover].push_back(Time);
m_PreTreatedData->SetCoreT(clover,crystal,Time);
@@ -221,7 +224,7 @@ void TGeTAMUPhysics::PreTreat(){
clover = m_EventData->GetSegmentCloverNbrE(i);
segment = m_EventData->GetSegmentSegmentNbrE(i);
name = "GETAMU/D"+ NPL::itoa(clover)+"_SEG"+ NPL::itoa(segment);
- Energy = cal->ApplyCalibration(name+"_E", Eraw);
+ Energy = cal->ApplyCalibration(name+"_E", Eraw+Random->Rndm());
Singles_CloverMap_SegEN[clover].push_back(segment);
Singles_CloverMap_SegE[clover].push_back(Energy);
m_PreTreatedData->SetSegmentE(clover,segment,Energy);
@@ -235,7 +238,7 @@ void TGeTAMUPhysics::PreTreat(){
clover = m_EventData->GetSegmentCloverNbrT(i);
segment = m_EventData->GetSegmentSegmentNbrT(i);
name = "GETAMU/D"+ NPL::itoa(clover)+"_SEG"+ NPL::itoa(segment);
- Time = cal->ApplyCalibration(name+"_T", Traw);
+ Time = cal->ApplyCalibration(name+"_T", Traw+Random->Rndm());
Singles_CloverMap_CryTN[clover].push_back(segment);
Singles_CloverMap_CryT[clover].push_back(Time);
m_PreTreatedData->SetSegmentT(clover,segment,Time);