ndmspc  0.20250128.0
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 
10 ClassImp(Ndmspc::Ndh::HnSparse);
12 
13 namespace Ndmspc {
14 namespace Ndh {
15 
16 HnSparse::HnSparse() : THnSparse()
17 {
21 }
22 
23 HnSparse::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 
32 Bool_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 
105 bool 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 
145 void 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