CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
LA_Fitter.cc
Go to the documentation of this file.
3 
4 #include <cmath>
5 #include <boost/foreach.hpp>
6 #include <boost/algorithm/string/erase.hpp>
7 #include <TF1.h>
8 
11  for(Book::iterator it = book.begin(".*"+method(WIDTH)); it!=book.end(); ++it) {
12  it->second->SetTitle("Mean Cluster Width;tan#theta_{t}");
13  TH1* const p = it->second;
14  if(p->GetEntries() < 400) { delete p; book[it->first]=0; continue;}
15  p->SetTitle(";tan#theta_{t};");
16  const float min = p->GetMinimum();
17  const float max = p->GetMaximum();
18  float xofmin = p->GetBinCenter(p->GetMinimumBin()); if( xofmin>0.0 || xofmin<-0.15) xofmin = -0.05;
19  const float xofmax = p->GetBinCenter(p->GetMaximumBin());
20 
21  TF1* const fit = new TF1("LA_profile_fit","[2]*(TMath::Abs(x-[0]))+[1]",-1,1);
22  fit->SetParLimits(0,-0.15,0.01);
23  fit->SetParLimits(1, 0.6*min, 1.25*min );
24  fit->SetParLimits(2,0.1,10);
25  fit->SetParameters( xofmin, min, (max-min) / fabs( xofmax - xofmin ) );
26 
27  int badfit = p->Fit(fit,"IEQ","",-.5,.3);
28  if( badfit ) badfit = p->Fit(fit,"IEQ","", -.46,.26);
29  if( badfit ) { book.erase(it); }
30  }
31 }
32 
35  for(Book::iterator it = book.begin(".*_all"); it!=book.end(); ++it) {
36  const std::string base = boost::erase_all_copy(it->first,"_all");
37 
38  std::vector<Book::iterator> rebin_hists;
39  Book::iterator all = it; rebin_hists.push_back(all);
40  Book::iterator w1 = book.find(base+"_w1"); rebin_hists.push_back(w1);
41  Book::iterator var_w2 = book.find(base+method(AVGV2,0)); rebin_hists.push_back(var_w2);
42  Book::iterator var_w3 = book.find(base+method(AVGV3,0)); rebin_hists.push_back(var_w3);
43 
44  const unsigned rebin = std::max( var_w2==book.end() ? 0 : find_rebin(var_w2->second),
45  var_w3==book.end() ? 0 : find_rebin(var_w3->second) );
46  BOOST_FOREACH(Book::iterator it, rebin_hists) if(it!=book.end()) it->second->Rebin( rebin>1 ? rebin<7 ? rebin : 6 : 1);
47 
48  TH1* const prob_w1 = w1==book.end() ? 0 : subset_probability( base+method(PROB1,0) ,w1->second,all->second);
49  TH1* const rmsv_w2 = var_w2==book.end() ? 0 : rms_profile( base+method(RMSV2,0), (TProfile*const)var_w2->second);
50  TH1* const rmsv_w3 = var_w3==book.end() ? 0 : rms_profile( base+method(RMSV3,0), (TProfile*const)var_w3->second);
51 
52  std::vector<TH1*> fit_hists;
53  if(prob_w1) {
54  book.book(base+method(PROB1,0),prob_w1);
55  fit_hists.push_back(prob_w1); prob_w1->SetTitle("Width==1 Probability;tan#theta_{t}-(dx/dz)_{reco}");
56  }
57  if(var_w2!=book.end()) {
58  book.book(base+method(RMSV2,0),rmsv_w2);
59  fit_hists.push_back(var_w2->second); var_w2->second->SetTitle("Width==2 Mean Variance;tan#theta_{t}-(dx/dz)_{reco}");
60  fit_hists.push_back(rmsv_w2); rmsv_w2->SetTitle("Width==2 RMS Variance;tan#theta_{t}-(dx/dz)_{reco}");
61  }
62  if(var_w3!=book.end()) {
63  book.book(base+method(RMSV3,0),rmsv_w3);
64  fit_hists.push_back(var_w3->second); var_w3->second->SetTitle("Width==3 Mean Variance;tan#theta_{t}-(dx/dz)_{reco}");
65  fit_hists.push_back(rmsv_w3); rmsv_w3->SetTitle("Width==3 RMS Variance;tan#theta_{t}-(dx/dz)_{reco}");
66  }
67 
68  if(!fit_hists.size()) continue;
69  const unsigned bins = fit_hists[0]->GetNbinsX();
70  const unsigned guess = fit_hists[0]->FindBin(0);
71  const std::pair<unsigned,unsigned> range(guess-bins/30,guess+bins/30);
72 
73  BOOST_FOREACH(TH1*const hist, fit_hists) {
74  TH1*const chi2 = SymmetryFit::symmetryChi2(hist,range);
75  if(chi2) {book.book(chi2->GetName(),chi2); chi2->SetTitle("Symmetry #chi^{2};tan#theta_{t}-(dx/dz)_{reco}");}
76  }
77  }
78 }
79 
80 unsigned LA_Filler_Fitter::
81 find_rebin(const TH1* const hist) {
82  const double mean = hist->GetMean();
83  const double rms = hist->GetRMS();
84  const int begin = std::min( 1, hist->GetXaxis()->FindFixBin(mean-rms));
85  const int end = std::max(hist->GetNbinsX(), hist->GetXaxis()->FindFixBin(mean+rms)) + 1;
86  unsigned current_hole(0), max_hole(0);
87  for(int i=begin; i<end; i++) {
88  if(!hist->GetBinError(i)) current_hole++;
89  else if(current_hole) {max_hole = std::max(current_hole,max_hole); current_hole=0;}
90  }
91  return max_hole+1;
92 }
93 
95 rms_profile(const std::string name, const TProfile* const prof) {
96  const int bins = prof->GetNbinsX();
97  TH1* const rms = new TH1F(name.c_str(),"",bins, prof->GetBinLowEdge(1), prof->GetBinLowEdge(bins) + prof->GetBinWidth(bins) );
98  for(int i = 1; i<=bins; i++) {
99  const double Me = prof->GetBinError(i);
100  const double neff = prof->GetBinEntries(i); //Should be prof->GetBinEffectiveEntries(i);, not availible this version ROOT. This is only ok for unweighted fills
101  rms->SetBinContent(i, Me*sqrt(neff) );
102  rms->SetBinError(i, Me/sqrt(2) );
103  }
104  return rms;
105 }
106 
108 subset_probability(const std::string name, const TH1* const subset, const TH1* const total) {
109  const int bins = subset->GetNbinsX();
110  TH1* const prob = new TH1F(name.c_str(),"",bins, subset->GetBinLowEdge(1), subset->GetBinLowEdge(bins) + subset->GetBinWidth(bins) );
111  for(int i = 1; i<=bins; i++) {
112  const double s = subset->GetBinContent(i);
113  const double T = total->GetBinContent(i);
114  const double B = T-s;
115 
116  const double p = T? s/T : 0;
117  const double perr = T? ( (s&&B)? sqrt(s*s*B+B*B*s)/(T*T) : 1/T ) : 0;
118 
119  prob->SetBinContent(i,p);
120  prob->SetBinError(i,perr);
121  }
122  return prob;
123 }
tuple base
Main Program
Definition: newFWLiteAna.py:92
int i
Definition: DBlmapReader.cc:9
static void fit(Book &book)
double_binary B
Definition: DDStreamer.cc:234
static TH1 * subset_probability(const std::string name, const TH1 *const , const TH1 *const )
Definition: LA_Fitter.cc:108
static void fit_width_profile(Book &)
Definition: LA_Fitter.cc:10
static TH1 * symmetryChi2(std::string, const std::vector< TH1 * > &, const std::pair< unsigned, unsigned >)
Definition: SymmetryFit.cc:6
const T & max(const T &a, const T &b)
T sqrt(T t)
Definition: SSEVec.h:48
iterator begin(string_t re=".*")
Definition: Book.h:48
#define end
Definition: vmac.h:37
static std::string method(Method m, bool fit=true)
static unsigned find_rebin(const TH1 *const)
Definition: LA_Fitter.cc:81
static TH1 * rms_profile(const std::string, const TProfile *const)
Definition: LA_Fitter.cc:95
static void make_and_fit_symmchi2(Book &)
Definition: LA_Fitter.cc:34
#define begin
Definition: vmac.h:30
Definition: Book.h:16
long double T
iterator find(string_t name, string_t re=".*")
Definition: Book.h:52
void erase(string_t name)
Definition: Book.h:57
boost::filter_iterator< match_name, book_t::iterator > iterator
Definition: Book.h:46
iterator end(string_t re=".*")
Definition: Book.h:50
TH1 * book(string_t name, TH1 *const hist)
Definition: Book.h:42