File indexing completed on 2023-03-17 10:39:26
0001
0002
0003
0004
0005 #include "PlotMilleMonitor.h"
0006 #include <TH1.h>
0007
0008
0009 #include <TF1.h>
0010 #include <TObjArray.h>
0011
0012 #include <TError.h>
0013 #include <TFile.h>
0014 #include <TROOT.h>
0015
0016
0017
0018
0019 #include "GFUtils/GFHistManager.h"
0020
0021
0022 PlotMilleMonitor::PlotMilleMonitor(const char *fileLegendList)
0023 : fHistManager(new GFHistManager)
0024 {
0025
0026
0027 if (!this->OpenFilesLegends(fileLegendList)) {
0028 ::Error("PlotMilleMonitor", "Problem opening files from '%s'", fileLegendList);
0029 }
0030 fHistManager->SameWithStats(true);
0031 fHistManager->SetLegendX1Y1X2Y2(0.14, 0.7, 0.45, 0.9);
0032
0033 }
0034
0035
0036 PlotMilleMonitor::~PlotMilleMonitor()
0037 {
0038 delete fHistManager;
0039
0040 for (unsigned int i = 0; i < fFileLegends.size(); ++i) delete fFileLegends[i].first;
0041
0042 }
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059
0060
0061
0062
0063
0064
0065
0066 void PlotMilleMonitor::DrawAllByHit(const char *xy, Option_t * option)
0067 {
0068 bool wasBatch = fHistManager->SetBatch();
0069
0070 TString opt(option);
0071 this->DrawResidualsByHit("resid", xy, opt);
0072 opt += " add";
0073 this->DrawResidualsByHit("reduResid", xy, opt);
0074 this->DrawResidualsByHit("sigma", xy, opt);
0075 this->DrawResidualsByHit("angle", xy, opt);
0076
0077 fHistManager->SetBatch(wasBatch);
0078 if (!wasBatch) fHistManager->Draw();
0079 }
0080
0081
0082
0083
0084
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097
0098
0099
0100
0101
0102
0103
0104
0105
0106
0107
0108
0109
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127
0128 void PlotMilleMonitor::DrawResidualsByHit(const char *histName, const char *xy, Option_t *option)
0129 {
0130 Int_t layer = this->PrepareAdd(TString(option).Contains("add", TString::kIgnoreCase));
0131
0132 for (unsigned int iFile = 0; iFile < fFileLegends.size(); ++iFile) {
0133 this->AddResidualsByHit(histName, fFileLegends[iFile], layer, xy, option);
0134 }
0135
0136 fHistManager->Draw();
0137 }
0138
0139
0140
0141 Int_t PlotMilleMonitor::AddResidualsByHit(const char *histName,
0142 std::pair<TFile*,TString> &fileLeg, Int_t layer,
0143 const char *xy, Option_t *option)
0144 {
0145
0146
0147
0148
0149 const TString opt(option);
0150 const bool addSumary = opt.Contains("sum", TString::kIgnoreCase);
0151 const bool onlySummary = opt.Contains("sumonly", TString::kIgnoreCase);
0152 const bool norm = opt.Contains("norm", TString::kIgnoreCase);
0153 const bool gaus = opt.Contains("gaus", TString::kIgnoreCase);
0154
0155 const unsigned int nHitsMax = 30;
0156 TH1 *hMean = 0, *hRms = 0;
0157 if (addSumary) {
0158 hMean = new TH1F(this->Unique(Form("%sMean", histName)),
0159 Form("mean of %s;N(hit)", histName), nHitsMax, 0, nHitsMax);
0160 hRms = new TH1F(this->Unique(Form("%sRms", histName)),
0161 Form("%s of %s;N(hit)", (gaus ? "#sigma" : "RMS"), histName),
0162 nHitsMax, 0, nHitsMax);
0163 }
0164 for (unsigned int i = 0; i < nHitsMax; ++i) {
0165 const TString completeName(Form("residuals/%s/%s_%u", xy, histName, i));
0166 TH1 *h = 0;
0167 fileLeg.first->GetObject(completeName, h);
0168 if (!h) {
0169 ::Warning("PlotMilleMonitor::DrawResidualsByHit", "No hist '%s'!", completeName.Data());
0170 } else {
0171 if (gaus) {
0172 h->Fit("gaus", "Q0L");
0173 h->GetFunction("gaus")->ResetBit(TF1::kNotDraw);
0174 }
0175 if (addSumary) {
0176 if (i == 0) {
0177 hMean->SetYTitle(Form("<%s>", h->GetXaxis()->GetTitle()));
0178 hRms->SetYTitle(Form("%s(%s)", (gaus ? "#sigma" : "RMS"), h->GetXaxis()->GetTitle()));
0179 }
0180 if (gaus) {
0181 TF1 *func = static_cast<TF1*>(h->GetFunction("gaus"));
0182 hMean->SetBinContent(i+1, func->GetParameter(1));
0183 hMean->SetBinError (i+1, func->GetParError(1));
0184 hRms-> SetBinContent(i+1, func->GetParameter(2));
0185 hRms-> SetBinError (i+1, func->GetParError(2));
0186 } else {
0187 hMean->SetBinContent(i+1, h->GetMean());
0188 hMean->SetBinError (i+1, h->GetMeanError());
0189 hRms-> SetBinContent(i+1, h->GetRMS());
0190 hRms-> SetBinError (i+1, h->GetRMSError());
0191 }
0192 }
0193 if (onlySummary) {
0194 delete h;
0195 } else {
0196 if (norm && h->GetEntries()) h->Scale(1./h->GetEntries());
0197 fHistManager->AddHistSame(h, layer, i, fileLeg.second);
0198 }
0199 }
0200 }
0201 if (hMean) fHistManager->AddHistSame(hMean, layer + (1 * !onlySummary), 0, fileLeg.second);
0202 if (hRms) fHistManager->AddHistSame(hRms, layer + (1 * !onlySummary), 1, fileLeg.second);
0203
0204 return ((addSumary && onlySummary) || (!addSumary) ? 1 : 2);
0205 }
0206
0207
0208 Int_t PlotMilleMonitor::PrepareAdd(bool addPlots)
0209 {
0210 if (addPlots) {
0211 return fHistManager->GetNumLayers();
0212 } else {
0213 fHistManager->Clear();
0214 return 0;
0215 }
0216 }
0217
0218
0219 TString PlotMilleMonitor::Unique(const char *name) const
0220 {
0221 if (!gROOT->FindObject(name)) return name;
0222
0223 UInt_t i = 1;
0224 while (gROOT->FindObject(Form("%s_%u", name, i))) ++i;
0225
0226 return Form("%s_%u", name, i);
0227 }
0228
0229
0230
0231 bool PlotMilleMonitor::OpenFilesLegends(const char *fileLegendList)
0232 {
0233 bool allOk = true;
0234
0235 TObjArray *fileLegPairs = TString(fileLegendList).Tokenize(",");
0236 for (Int_t iF = 0; iF < fileLegPairs->GetEntriesFast(); ++iF) {
0237 TObjArray *aFileLegPair = TString(fileLegPairs->At(iF)->GetName()).Tokenize("=");
0238
0239 const char *legend = "";
0240 if (aFileLegPair->GetEntriesFast() >= 2) {
0241 if (aFileLegPair->GetEntriesFast() > 2) {
0242 ::Error("TifResidualOverlay::OpenFilesLegends",
0243 "File-legend pair %s: %d (>2) '=' separated parts!",
0244 fileLegPairs->At(iF)->GetName(), aFileLegPair->GetEntriesFast());
0245 }
0246 legend = aFileLegPair->At(1)->GetName();
0247 } else if (aFileLegPair->GetEntriesFast() < 1) {
0248 continue;
0249 }
0250
0251 TFile *file = TFile::Open(aFileLegPair->At(0)->GetName());
0252 if (!file) {
0253 allOk = false;
0254 } else {
0255
0256 fFileLegends.push_back(std::make_pair(file, legend));
0257 }
0258 delete aFileLegPair;
0259 }
0260 delete fileLegPairs;
0261
0262 return allOk;
0263 }