ThresholdFitAlgorithm2.cpp

00001 #include "ThresholdFitAlgorithm2.h"
00002 #include "FitAlgorithmMap.h"
00003 #include "FitterMode.h"
00004 #include "SctData/mVThresholdVariable.h"
00005 #include "SctData/mVfromTrimTargetThresholdVariable.h"
00006 #include "SctData/ErfcFitObject.h"
00007 #include "SctData/StandardDefects.h"
00008 #include <TH1.h>
00009 
00010 using namespace SctData;
00011 
00012 namespace SctFitter {
00013     
00014     
00015 ThresholdFitAlgorithm2::ThresholdFitAlgorithm2() throw() {
00016     FitterMode mode = getMode();
00017     mode.fitNone(); 
00018     mode.fitChips(); 
00019     mode.fitChannels() ; 
00020     mode.doSummary(false);     
00021     setMode(mode);
00022 }
00023 
00024 bool ThresholdFitAlgorithm2::putInMap() throw() {    
00025     return FitAlgorithmMap::instance().setAlgorithm("ThresholdAnalytic", auto_ptr<FitAlgorithm>(new ThresholdFitAlgorithm2));
00026 }
00027     
00028 bool ThresholdFitAlgorithm2::inMap = putInMap();
00029 
00030 void ThresholdFitAlgorithm2::doFit(const TH1& hist, FitObject& fitObject,
00031                    const ModuleElement& element, DefectList& defects) const throw (LogicError) {
00032     
00033     double amp = hist.GetBinContent(1);
00034     unsigned int ampCount = 1;
00035     double mean = 0;
00036     double meanx2 = 0;
00037     double norm = 0;
00038     double last = hist.GetBinContent(1);
00039     bool flatTop = true;
00040     
00041     for (unsigned int i=2; i<=hist.GetNbinsX(); ++i) {
00042     double bin = hist.GetBinContent(i);
00043     double bincentre = hist.GetBinLowEdge(i);
00044     double error = hist.GetBinError(i);
00045     double diff = bin-last;
00046     norm += diff;
00047     if (flatTop && error < (amp/ampCount - bin)) {
00048         //cout << "diff>error: " << amp << " error: " << error << " last: "<< last << " bin: " << i << " bincentre" << bincentre << endl; 
00049         flatTop = false;
00050     }
00051     if (flatTop) {
00052         amp += bin;
00053         ++ampCount;
00054     }
00055     mean += diff * bincentre;
00056     meanx2 += diff*bincentre* bincentre;
00057     last = bin;
00058     }
00059  
00060     //cout << "amp: " << amp << " ampCount: " << ampCount << " mean: " << mean << " meanx2 " << meanx2 << endl;
00061     
00062     amp /= ampCount;
00063     //mean /= 1-hist.GetNbinsX();
00064     mean /= norm;
00065     //double sigma = sqrt(meanx2 /(1-hist.GetNbinsX()) - mean * mean);
00066     double sigma = sqrt(meanx2/norm - mean * mean);
00067     
00068     //cout << "amp: " << amp << " mean: " << mean << " sigma: " << sigma << endl;;
00069     
00070     
00071     fitObject.setParameter(0, amp);
00072     fitObject.setParameter(1, mean);
00073     fitObject.setParameter(2, sigma);
00074     fitObject.setParError(0, 0);
00075     fitObject.setParError(1, 0);
00076     fitObject.setParError(2, 0);    
00077 }
00078 
00079     
00080 auto_ptr<FitObject> ThresholdFitAlgorithm2::getPrototype() const throw() {
00081     return auto_ptr<FitObject> (new ErfcFitObject());
00082 }
00083     
00084 void ThresholdFitAlgorithm2::guessParameters(const TH1& hist, FitObject& fitOb) const throw (LogicError, MathsError) {
00085 }
00086 
00087 
00088 void ThresholdFitAlgorithm2::checkForDefects(const FitObject& fo, const ModuleElement& element, 
00089                          DefectList& defects) const {
00090   if (fo.getNDF() &&  fo.getChiSquared()/fo.getNDF() > StandardDefects::BADFIT.getParameter() ) {
00091     cout << "Chi-squared: " << fo.getChiSquared() << " NDF: " << fo.getNDF() << endl;
00092       defects.addDefect(Defect(StandardDefects::BADFIT,element));
00093   }
00094      
00095 }
00096 
00097 void ThresholdFitAlgorithm2::checkForDefects(const TH1& hist, const ModuleElement& element, 
00098                          DefectList& defects) const throw (LogicError) {
00099 }
00100 
00101 void ThresholdFitAlgorithm2::createSummaryHistograms(FitScanResult& fits) const throw() {
00102 }
00103 }

Generated on Mon Feb 6 14:01:35 2006 for SCT DAQ/DCS Software - C++ by  doxygen 1.4.6