https://github.com/root-project/root
Tip revision: 729005709a50b1e2e8b5400a8a20c209af5c7cbe authored by Axel Naumann on 29 September 2022, 13:03:06 UTC
"Update ROOT version files to v6.24/08."
"Update ROOT version files to v6.24/08."
Tip revision: 7290057
dt_MakeRef.C
#include "dt_RunDrawTest.C"
#include "TClassTable.h"
#include "TDirectory.h"
#include "TFile.h"
#include "TH1.h"
#include "TH2.h"
#include "TSystem.h"
#ifndef __CINT__
#include "Event.h"
#endif
R__LOAD_LIBRARY(libEvent)
template <class Hist>
Hist *RefClone(Hist* orig) {
Hist *cloned = (Hist*)orig->Clone();
TString name = orig->GetName();
name.Prepend("ref");
cloned->SetName(name);
cloned->Reset();
return cloned;
};
template <class Hist>
Hist* RefClone(TDirectory* from, const char* name) {
Hist * orig = dynamic_cast<Hist*>(from->Get(name));
if (!orig) {
TObject *obj = from->Get(name);
if (obj) {
cerr << "Type of " << name << " is incorrect. It is "
<< obj->IsA()->GetName() << " instead of " << Hist::Class()->GetName() << endl;
} else {
cerr << "Missing " << name << " from " << from->GetName() << endl;
}
return 0;
}
return RefClone(orig);
}
void MakeHisto(TTree *tree, TDirectory* To) {
cout << "Generating histograms from TTree::Draw" << endl;
TDirectory* where = GenerateDrawHist(tree);
To->cd();
Event *event = new Event();
tree->SetBranchAddress("event",&event);
//We make clones of the generated histograms
//We set new names and reset the clones.
//We want to have identical histogram limits
TH1F *refNtrack = RefClone<TH1F>(where,"hNtrack");
TH1F *refGetNtrack = RefClone<TH1F>(where,"hGetNtrack");
TH1F *refNseg = RefClone<TH1F>(where,"hNseg");
TH1F *refTemp = RefClone<TH1F>(where,"hTemp");
TH1F *refHmean = RefClone<TH1F>(where,"hHmean");
TH1F *refHAxisMax = RefClone<TH1F>(where,"hHAxisMax");
TH1F *refHAxisGetMax = RefClone<TH1F>(where,"hHAxisGetMax");
TH1F *refHGetAxisGetMax = RefClone<TH1F>(where,"hHGetAxisGetMax");
TH1F *refHGetAxisMax = RefClone<TH1F>(where,"hHGetAxisMax");
TH1F *refGetHGetAxisMax = RefClone<TH1F>(where,"hGetHGetAxisMax");
TH1F *refGetRefHGetAxisMax = RefClone<TH1F>(where,"hGetRefHGetAxisMax");
TH1F *refPx = RefClone<TH1F>(where,"hPx");
TH1F *refPy = RefClone<TH1F>(where,"hPy");
TH1F *refPz = RefClone<TH1F>(where,"hPz");
TH1F *refRandom = RefClone<TH1F>(where,"hRandom");
TH1F *refMass2 = RefClone<TH1F>(where,"hMass2");
TH1F *refBx = RefClone<TH1F>(where,"hBx");
TH1F *refBy = RefClone<TH1F>(where,"hBy");
TH1F *refXfirst = RefClone<TH1F>(where,"hXfirst");
TH1F *refYfirst = RefClone<TH1F>(where,"hYfirst");
TH1F *refZfirst = RefClone<TH1F>(where,"hZfirst");
TH1F *refXlast = RefClone<TH1F>(where,"hXlast");
TH1F *refYlast = RefClone<TH1F>(where,"hYlast");
TH1F *refZlast = RefClone<TH1F>(where,"hZlast");
TH1F *refCharge = RefClone<TH1F>(where,"hCharge");
TH1F *refNpoint = RefClone<TH1F>(where,"hNpoint");
TH1F *refValid = RefClone<TH1F>(where,"hValid");
TH1F *refPointValue = RefClone<TH1F>(where,"hPointValue");
TH1F *refAlias = RefClone<TH1F>(where,"hAlias");
TH1F *refAliasSymbol = RefClone<TH1F>(where,"hAliasSymbol");
TH1F *refAliasSymbolFunc = RefClone<TH1F>(where,"hAliasSymbolFunc");
TH1F *refBool = RefClone<TH1F>(where,"hBool");
TH1F *refFullMatrix = RefClone<TH1F>(where,"hFullMatrix");
TH1F *refColMatrix = RefClone<TH1F>(where,"hColMatrix");
TH1F *refRowMatrix = RefClone<TH1F>(where,"hRowMatrix");
TH1F *refCellMatrix = RefClone<TH1F>(where,"hCellMatrix");
TH1F *refFullOper = RefClone<TH1F>(where,"hFullOper");
TH1F *refCellOper = RefClone<TH1F>(where,"hCellOper");
TH1F *refColOper = RefClone<TH1F>(where,"hColOper");
TH1F *refRowOper = RefClone<TH1F>(where,"hRowOper");
TH1F *refMatchRowOper = RefClone<TH1F>(where,"hMatchRowOper");
TH1F *refMatchColOper = RefClone<TH1F>(where,"hMatchColOper");
TH1F *refRowMatOper = RefClone<TH1F>(where,"hRowMatOper");
TH1F *refMatchDiffOper= RefClone<TH1F>(where,"hMatchDiffOper");
TH1F *refFullOper2 = RefClone<TH1F>(where,"hFullOper2");
TH1F *refClosestDistance = RefClone<TH1F>(where,"hClosestDistance");
TH1F *refClosestDistance2 = RefClone<TH1F>(where,"hClosestDistance2");
TH1F *refClosestDistance9 = RefClone<TH1F>(where,"hClosestDistance9");
TH1F *refClosestDistanceIndex = RefClone<TH1F>(where, "hClosestDistanceIndex");
TH2F *refPxInd = RefClone<TH2F>(where,"hPxInd");
TH1F *refSqrtNtrack = RefClone<TH1F>(where,"hSqrtNtrack");
TH1F *refShiftValid = RefClone<TH1F>(where,"hShiftValid");
TH1F *refAndValid = RefClone<TH1F>(where,"hAndValid");
TH1F *refString = RefClone<TH1F>(where,"hString");
TH1F *refStringSpace = RefClone<TH1F>(where,"hStringSpace");
TH1F *refAliasStr = RefClone<TH1F>(where,"hAliasStr");
TH1F *refPxBx = RefClone<TH1F>(where,"hPxBx");
TH1F *refPxBxWeight = RefClone<TH1F>(where,"hPxBxWeight");
TH1F *refTriggerBits = RefClone<TH1F>(where,"hTriggerBits");
TH1F *refTriggerBitsFunc = RefClone<TH1F>(where,"hTriggerBitsFunc");
TH1F *refFiltTriggerBits = RefClone<TH1F>(where,"hFiltTriggerBits");
TH1F *refTrackTrigger = RefClone<TH1F>(where,"hTrackTrigger");
TH1F *refFiltTrackTrigger = RefClone<TH1F>(where,"hFiltTrackTrigger");
TH1F *refBreit = RefClone<TH1F>(where,"hBreit");
TH1F *refAlt = RefClone<TH1F>(where,"hAlt");
TH1F *refSize = RefClone<TH1F>(where,"hSize");
TH1F *refSize2 = RefClone<TH1F>(where,"hSize2");
TH1F *refSumPx = RefClone<TH1F>(where,"hSumPx");
TH2F *refMaxPx = RefClone<TH2F>(where,"hMaxPx");
TH2F *refMinPx = RefClone<TH2F>(where,"hMinPx");
// Loop with user code on all events and fill the ref histograms
// The code below should produce identical results to the tree->Draw above
cout << "Recalculating the histograms with custom loop." << endl;
TClonesArray *tracks = event->GetTracks();
Int_t nev = (Int_t)tree->GetEntries();
Int_t i, ntracks, evmod,i0,i1, Nvertex;
Track *t;
EventHeader *head;
Int_t nbin = 0;
for (Int_t ev=0;ev<nev;ev++) {
nbin += tree->GetEntry(ev);
head = event->GetHeader();
evmod = head->GetEvtNum()%10;
refNtrack->Fill(event->GetNtrack());
refGetNtrack->Fill(event->GetNtrack());
refNseg->Fill(event->GetNseg());
refTemp->Fill(event->GetTemperature());
refHmean->Fill(event->GetHistogram()->GetMean());
refHAxisMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refHAxisGetMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refHGetAxisGetMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refHGetAxisMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refGetHGetAxisMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refGetRefHGetAxisMax->Fill(event->GetHistogram()->GetXaxis()->GetXmax());
refSqrtNtrack->Fill(sqrt(event->GetNtrack()));
if (!strcmp("type1",event->GetType()))
refString->Fill(event->GetHeader()->GetEvtNum());
if (strstr(event->GetType(),"1")) {
refString->Fill(event->GetHeader()->GetEvtNum());
}
if (strcmp("Event Histogram",event->GetHistogram()->GetTitle())==0) {
refStringSpace->Fill(1);
}
refAliasStr->Fill(strstr(event->GetType(),"1")!=0);
refBool->Fill(event->IsValid());
Nvertex = event->GetNvertex();
for(i0=0;i0<Nvertex;i0++) {
refClosestDistance->Fill(event->GetClosestDistance(i0));
}
if (Nvertex>2) refClosestDistance2->Fill(event->GetClosestDistance(2));
if (Nvertex>9) refClosestDistance9->Fill(event->GetClosestDistance(9));
refClosestDistanceIndex->Fill(event->GetClosestDistance(Nvertex/2));
for(i0=0;i0<4;i0++) {
for(i1=0;i1<4;i1++) {
refFullMatrix->Fill(event->GetMatrix(i0,i1));
int i2 = i0*4+i1;
refAlt->Fill( event->GetMatrix(i0,i1) - ( (i2<Nvertex) ? event->GetClosestDistance(i2) : 0 ) );
}
refColMatrix->Fill(event->GetMatrix(i0,0));
refRowMatrix->Fill(event->GetMatrix(1,i0)); // done here because the matrix is square!
}
refCellMatrix->Fill(event->GetMatrix(2,2));
TBits bits = event->GetTriggerBits();
Int_t nbits = bits.GetNbits();
Int_t ncx = refTriggerBits->GetXaxis()->GetNbins();
Int_t nextbit = -1;
while(1) {
nextbit = bits.FirstSetBit(nextbit+1);
if (nextbit >= nbits) break;
if (nextbit > ncx) refTriggerBits->Fill(ncx+1);
else refTriggerBits->Fill(nextbit);
if (nextbit > ncx) refTriggerBitsFunc->Fill(ncx+1);
else refTriggerBitsFunc->Fill(nextbit);
}
if (bits.TestBitNumber(28)) refFiltTriggerBits->Fill(nbits);
ntracks = event->GetNtrack();
refSize->Fill(ntracks);
refSize->Fill(ntracks);
refSize2->Fill(ntracks);
refSize2->Fill(ntracks);
if ( 5 < ntracks ) {
t = (Track*)tracks->UncheckedAt(5);
for(i0=0;i0<4;i0++) {
for(i1=0;i1<4;i1++) {
}
refColOper->Fill( event->GetMatrix(i0,1) - t->GetVertex(1) );
refRowOper->Fill( event->GetMatrix(2,i0) - t->GetVertex(2) );
}
for(i0=0;i0<3;i0++) {
refMatchRowOper->Fill( event->GetMatrix(2,i0) - t->GetVertex(i0) );
refMatchDiffOper->Fill( event->GetMatrix(i0,2) - t->GetVertex(i0) );
}
refCellOper->Fill( event->GetMatrix(2,1) - t->GetVertex(1) );
}
Double_t sumPx = 0;
Double_t maxPx = 0.0;
Double_t maxPy = 0.0;
Double_t minPx = 5.0;
Double_t minPy = 5.0;
for (i=0;i<ntracks;i++) {
t = (Track*)tracks->UncheckedAt(i);
sumPx += t->GetPx();
if (t->GetPy() > 1.0) {
if (t->GetPx() > maxPx) maxPx = t->GetPx();
if (t->GetPx() < minPx) minPx = t->GetPx();
}
if (t->GetPy() > maxPy) maxPy = t->GetPy();
if (t->GetPy() < minPy) minPy = t->GetPy();
if (evmod == 0) refPx->Fill(t->GetPx());
if (evmod == 0) refPy->Fill(t->GetPy());
if (evmod == 0) refPz->Fill(t->GetPz());
if (evmod == 1) refRandom->Fill(t->GetRandom(),3);
if (evmod == 1) refMass2->Fill(t->GetMass2());
if (evmod == 1) refBx->Fill(t->GetBx());
if (evmod == 1) refBy->Fill(t->GetBy());
if (evmod == 2) refXfirst->Fill(t->GetXfirst());
if (evmod == 2) refYfirst->Fill(t->GetYfirst());
if (evmod == 2) refZfirst->Fill(t->GetZfirst());
if (evmod == 3) refXlast->Fill(t->GetXlast());
if (evmod == 3) refYlast->Fill(t->GetYlast());
if (evmod == 3) refZlast->Fill(t->GetZlast());
if (t->GetPx() < 0) {
refCharge->Fill(t->GetCharge());
refNpoint->Fill(t->GetNpoint());
refValid->Fill(t->GetValid());
}
Int_t valid = t->GetValid();
refShiftValid->Fill(valid << 4);
refShiftValid->Fill( (valid << 4) >> 2 );
if (event->GetNvertex()>10 && event->GetNseg()<=6000) {
refAndValid->Fill( t->GetValid() & 0x1 );
}
Track * t2 = (Track*)tracks->At(t->GetNpoint()/6);
if (t2 && t2->GetPy()>0) {
refPxInd->Fill(t2->GetPy(),t->GetPx());
}
float Bx,By;
Bx = t->GetBx();
By = t->GetBy();
if ((Bx>.15) || (By<=-.15)) refPxBx->Fill(t->GetPx());
double weight = Bx*Bx*(Bx>.15) + By*By*(By<=-.15);
if (weight) refPxBxWeight->Fill(t->GetPx(),weight);
if (i<4) {
for(i1=0;i1<3;i1++) { // 3 is the min of the 2nd dim of Matrix and Vertex
refFullOper ->Fill( event->GetMatrix(i,i1) - t->GetVertex(i1) );
refFullOper2->Fill( event->GetMatrix(i,i1) - t->GetVertex(i1) );
refRowMatOper->Fill( event->GetMatrix(i,2) - t->GetVertex(i1) );
}
refMatchColOper->Fill( event->GetMatrix(i,2) - t->GetVertex(1) );
}
for(i1=0; i1<t->GetN(); i1++) {
refPointValue->Fill( t->GetPointValue(i1) );
}
TBits bits = t->GetTriggerBits();
Int_t nbits = bits.GetNbits();
Int_t ncx = refTrackTrigger->GetXaxis()->GetNbins();
Int_t nextbit = -1;
while(1) {
nextbit = bits.FirstSetBit(nextbit+1);
if (nextbit >= nbits) break;
if (nextbit > ncx) refTrackTrigger->Fill(ncx+1);
else refTrackTrigger->Fill(nextbit);
}
if (bits.TestBitNumber(5)) refFiltTrackTrigger->Fill(t->GetPx());
refBreit->Fill(TMath::BreitWigner(t->GetPx(),3,2));
refAlias->Fill(head->GetEvtNum()*6+t->GetPx()*t->GetPy());
refAliasSymbol->Fill(t->GetPx()+t->GetPy());
refAliasSymbolFunc->Fill(t->GetPx()+t->GetPy());
}
refSumPx->Fill(sumPx);
if (maxPx > 0) {
refMaxPx->Fill(maxPy,maxPx);
}
if (minPx < 5.0) {
refMinPx->Fill(minPy,minPx);
}
}
delete event;
Event::Reset();
}
void dt_MakeRef(const char* from, Int_t verboseLevel = 2) {
SetVerboseLevel(verboseLevel);
if (!TClassTable::GetDict("Event")) {
gSystem->Load("libEvent");
gHasLibrary = kTRUE;
} else {
gHasLibrary = kTRUE;
}
gROOT->GetList()->Delete();
TFile *hfile = new TFile(from);
TTree *tree = (TTree*)hfile->Get("T");
TFile* f= new TFile("dt_reference.root","recreate");
MakeHisto(tree,f);
f->Write();
delete f;
delete hfile;
gROOT->cd();
cout << "Checking histograms" << endl;
Compare(gDirectory);
// gROOT->GetList()->Delete();
}