FairRoot
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
FairMQPixelFileSinkBin.cxx
Go to the documentation of this file.
1 /********************************************************************************
2  * Copyright (C) 2014 GSI Helmholtzzentrum fuer Schwerionenforschung GmbH *
3  * *
4  * This software is distributed under the terms of the *
5  * GNU Lesser General Public Licence (LGPL) version 3, *
6  * copied verbatim in the file "LICENSE" *
7  ********************************************************************************/
15 #include "FairMQPixelFileSinkBin.h"
16 
17 #include "FairEventHeader.h"
18 #include "PixelHit.h"
19 #include "PixelPayload.h"
20 
21 #include <FairMQLogger.h>
22 #include <TClonesArray.h>
23 #include <TFile.h>
24 #include <TFolder.h>
25 #include <TList.h>
26 #include <TObjString.h>
27 #include <TObject.h>
28 #include <TTree.h>
29 #include <TVector3.h>
30 
31 using namespace std;
32 
34  : FairMQDevice()
35  , fInputChannelName("data-in")
36  , fAckChannelName("")
37  , fFileName()
38  , fTreeName()
39 
40  , fBranchNames()
41  , fClassNames()
42  , fFileOption()
43  , fFlowMode(false)
44  , fWrite(false)
45 
46  , fOutFile(nullptr)
47  , fTree(nullptr)
48  , fNObjects(0)
49  , fOutputObjects(new TObject*[1000])
50  , fFolder(nullptr)
51 {}
52 
54 {
55  fFileName = fConfig->GetValue<std::string>("file-name");
56  fClassNames = fConfig->GetValue<std::vector<std::string>>("class-name");
57  fBranchNames = fConfig->GetValue<std::vector<std::string>>("branch-name");
58  fInputChannelName = fConfig->GetValue<std::string>("in-channel");
59  fAckChannelName = fConfig->GetValue<std::string>("ack-channel");
60 
61  LOG(info) << "SHOULD CREATE THE FILE AND TREE";
62  fFileOption = "RECREATE";
63  fTreeName = "cbmsim";
64 
65  fOutFile = TFile::Open(fFileName.c_str(), fFileOption.c_str());
66 
67  fTree = new TTree(fTreeName.c_str(), "/cbmout");
68 
69  fFolder = new TFolder("cbmout", "Main Output Folder");
70  TFolder* foldEventHeader = fFolder->AddFolder("EvtHeader", "EvtHeader");
71  TFolder* foldPixel = fFolder->AddFolder("Pixel", "Pixel");
72 
73  TList branchNameList;
74 
75  for (fNObjects = 0; fNObjects < fBranchNames.size(); fNObjects++) {
76  if (fClassNames[fNObjects].find("TClonesArray(") == 0) {
77  fClassNames[fNObjects] = fClassNames[fNObjects].substr(13, fClassNames[fNObjects].length() - 12 - 2);
78  fOutputObjects[fNObjects] = new TClonesArray(fClassNames[fNObjects].c_str());
79  fTree->Branch(fBranchNames[fNObjects].c_str(), "TClonesArray", &fOutputObjects[fNObjects]);
80  foldPixel->Add(fOutputObjects[fNObjects]);
81  branchNameList.AddLast(new TObjString(fBranchNames[fNObjects].c_str()));
82  } else if (fClassNames[fNObjects].find("FairEventHeader") == 0) {
83  fOutputObjects[fNObjects] = new FairEventHeader();
84  fTree->Branch(fBranchNames[fNObjects].c_str(), "FairEventHeader", &fOutputObjects[fNObjects]);
85  foldEventHeader->Add(fOutputObjects[fNObjects]);
86  branchNameList.AddLast(new TObjString(fBranchNames[fNObjects].c_str()));
87  } else {
88  LOG(error) << "!!! Unknown output object \"" << fClassNames[fNObjects] << "\" !!!";
89  }
90  }
91 
92  fFolder->Write();
93  branchNameList.Write("BranchList", TObject::kSingleKey);
94  branchNameList.Delete();
95 
96  OnData(fInputChannelName, &FairMQPixelFileSinkBin::StoreData);
97 }
98 
99 bool FairMQPixelFileSinkBin::StoreData(FairMQParts& parts, int /*index*/)
100 {
101  if (parts.Size() == 0)
102  return true; // probably impossible, but still check
103 
104  // expecting even number of parts in the form: header,data,header,data,header,data and so on...
105  int nPPE = 2; // nof parts per event
106 
107  if (parts.Size() % nPPE >= 1)
108  LOG(info) << "received " << parts.Size() << " parts, will ignore last part!!!";
109 
110  for (int ievent = 0; ievent < parts.Size() / nPPE; ievent++) {
111  // the first part should be the event header
112  PixelPayload::EventHeader* payloadE =
113  static_cast<PixelPayload::EventHeader*>(parts.At(nPPE * ievent)->GetData());
114  // LOG(debug) << "GOT EVENT " << payloadE->fMCEntryNo << " OF RUN " << payloadE->fRunId << " (part " <<
115  // payloadE->fPartNo << ")";
116 
117  for (unsigned int ibr = 0; ibr < fBranchNames.size(); ibr++) {
118  if ("EventHeader." == fBranchNames[ibr]) {
119  ((FairEventHeader*)fOutputObjects[ibr])->SetRunId(payloadE->fRunId);
120  ((FairEventHeader*)fOutputObjects[ibr])->SetMCEntryNumber(payloadE->fMCEntryNo);
121  }
122  }
123 
124  // the second part should the TClonesArray with necessary data... now assuming Digi
125  PixelPayload::Hit* payloadH = static_cast<PixelPayload::Hit*>(parts.At(nPPE * ievent + 1)->GetData());
126  int hitArraySize = parts.At(nPPE * ievent + 1)->GetSize();
127  int nofHits = hitArraySize / sizeof(PixelPayload::Hit);
128 
129  for (unsigned int ibr = 0; ibr < fBranchNames.size(); ibr++) {
130  if ("PixelHits" == fBranchNames[ibr]) {
131  ((TClonesArray*)fOutputObjects[ibr])->Clear();
132  for (int ihit = 0; ihit < nofHits; ihit++) {
133  TVector3 pos(payloadH[ihit].posX, payloadH[ihit].posY, payloadH[ihit].posZ);
134  TVector3 posErr(payloadH[ihit].dposX, payloadH[ihit].dposY, payloadH[ihit].dposZ);
135  new ((*((TClonesArray*)fOutputObjects[ibr]))[ihit])
136  PixelHit(payloadH[ihit].fDetectorID, -1, pos, posErr);
137  // new ((*fHits)[fNHits]) PixelHit(detId,currentDigi->GetIndex(),pos,posErr);
138  }
139  }
140  }
141 
142  fTree->Fill();
143  }
144 
145  if (fAckChannelName != "") {
146  FairMQMessagePtr msg(NewMessage());
147  Send(msg, fAckChannelName);
148  }
149  return true;
150 }
151 
153 {
154  if (fTree) {
155  fTree->Write();
156  delete fTree;
157  }
158 
159  if (fOutFile) {
160  if (fOutFile->IsOpen()) {
161  fOutFile->Close();
162  }
163  delete fOutFile;
164  }
165 }
bool StoreData(FairMQParts &, int)