ndmspc 0.20250304.0
Loading...
Searching...
No Matches
HnSparse.cxx
1#include <TFile.h>
2#include <TTree.h>
3#include <TAxis.h>
4#include <TH2D.h>
5#include <THnSparse.h>
6
7#include "HnSparse.h"
8
10ClassImp(Ndmspc::Ndh::HnSparse);
12
13namespace Ndmspc {
14namespace Ndh {
15
16HnSparse::HnSparse() : THnSparse()
17{
21}
22
23HnSparse::HnSparse(const char * name, const char * title, Int_t dim, const Int_t * nbins, const Double_t * xmin,
24 const Double_t * xmax, Int_t chunksize)
25 : THnSparse(name, title, dim, nbins, xmin, xmax, chunksize)
26{
30}
31
32Bool_t HnSparse::Import(std::vector<Int_t> r, TString filename, TString objname, TString cacheDir)
33{
37
38 if (!cacheDir.IsNull()) TFile::SetCacheFileDir(cacheDir.Data(), 1, 1);
39
40 if (filename.IsNull()) {
41 Printf("Error: filename is empty !!!");
42 return kFALSE;
43 }
44 if (objname.IsNull()) {
45 Printf("Error: objname is empty !!!");
46 return kFALSE;
47 }
48
49 Printf("Opening file='%s' obj='%s' ...", filename.Data(), objname.Data());
50 TFile * f = TFile::Open(filename.Data());
51 if (f == nullptr) return kFALSE;
52
53 THnSparse * s = (THnSparse *)f->Get(objname.Data());
54
55 // s->Print();
56
57 TAxis * a;
58
59 TObjArray * newAxis = (TObjArray *)s->GetListOfAxes()->Clone();
60 for (Int_t iDim = 0; iDim < newAxis->GetEntries(); iDim++) {
61 a = (TAxis *)newAxis->At(iDim);
62 // Printf("%s %d %.2f %.2f", a->GetName(), a->GetNbins(), a->GetXmin(), a->GetXmax());
63 // a->Print();
64 if (std::find(r.begin(), r.end(), iDim) != r.end()) {
65 // Printf("%s %d %.2f %.2f", a->GetName(), a->GetNbins(), a->GetXmin(), a->GetXmax());
66 }
67 else {
68 /* v does not contain x */
69 // Printf("Reset %s %d %.2f %.2f", a->GetName(), a->GetNbins(), a->GetXmin(), a->GetXmax());
70 a->Set(1, a->GetXmin(), a->GetXmax());
71 }
72 }
73
74 Init(TString::Format("NDMSPC_%s", s->GetName()).Data(), "", newAxis, kTRUE);
75
76 Int_t dims[fNdimensions];
77 Int_t c[fNdimensions];
78 for (Int_t iDim = 0; iDim < newAxis->GetEntries(); iDim++) {
79 // a = s->GetAxis(iDim);
80 // Printf("%s %d %.2f %.2f", a->GetName(), a->GetNbins(), a->GetXmin(), a->GetXmax());
81 dims[iDim] = iDim;
82 c[iDim] = 1;
83 }
84 TFile * fileOut = TFile::Open(fOutputFileName.Data(), "RECREATE");
85 fTree = new TTree("ndh", "Ndh tree");
86 fTree->Branch("h", &s);
87
88 RecursiveLoop(s, 0, c, dims, r);
89
90 // Reset cuts
91 for (Int_t iDim = 0; iDim < GetNdimensions(); iDim++) {
92 GetAxis(iDim)->SetRange();
93 }
94
95 fTree->GetUserInfo()->Add(Clone());
96 fTree->Write();
97 fileOut->Close();
98
99 delete s;
100 f->Close();
101
102 return kTRUE;
103}
104
105bool HnSparse::RecursiveLoop(THnSparse * s, Int_t level, Int_t * coord, Int_t * dims, std::vector<Int_t> & r)
106{
110
111 if (level >= (Int_t)r.size()) return true;
112
113 // Printf("level=%d axis_id=%d", level, r[level]);
114
115 for (Int_t iBin = 1; iBin <= GetAxis(r[level])->GetNbins(); iBin++) {
116
117 coord[r[level]] = iBin;
118 s->GetAxis(r[level])->SetRange(iBin, iBin);
119 Bool_t finished = RecursiveLoop(s, level + 1, coord, dims, r);
120 if (finished) {
121 THnSparse * ss = (THnSparse *)s->ProjectionND(s->GetNdimensions(), dims, "O");
122
123 ss->SetName(GetName());
124 ss->SetEntries(1);
125 // GetBin(coord, kTRUE);
126 if (ss->GetNbins() > 0) {
127 SetBinContent(coord, ss->GetNbins());
128 Printf("level=%d axis_id=%d iBin=%d binsFilled=%lld", level, r[level], iBin, ss->GetNbins());
129 fTree->SetBranchAddress("h", &ss);
130 fTree->Fill();
131 }
132 else {
133 Printf("[NotFilled] level=%d axis_id=%d iBin=%d binsFilled=%lld", level, r[level], iBin, ss->GetNbins());
134 }
135 delete ss;
136 }
137 else {
138 Printf("level=%d axis_id=%d iBin=%d", level, r[level], iBin);
139 }
140 }
141
142 return false;
143}
144
145void HnSparse::ReserveBins(Long64_t nBins)
146{
150 Printf("Reserving %e bins ...", (Double_t)nBins);
151 Reserve(nBins);
152 Printf("%e bins reserved.", (Double_t)nBins);
153}
154} // namespace Ndh
155} // namespace Ndmspc
HnSparse object.
Definition HnSparse.h:19
TTree * fTree
Container.
Definition HnSparse.h:40
void ReserveBins(Long64_t nBins)
Definition HnSparse.cxx:145
Bool_t Import(std::vector< Int_t > r, TString filename, TString objname, TString cacheDir=gSystem->HomeDirectory())
Definition HnSparse.cxx:32
bool RecursiveLoop(THnSparse *s, Int_t level, Int_t *coord, Int_t *dims, std::vector< Int_t > &r)
Definition HnSparse.cxx:105
TString fOutputFileName
Output filename.
Definition HnSparse.h:41