1 | // %%%%%%%%%%%%% |
---|
2 | // gemc headers |
---|
3 | // %%%%%%%%%%%%% |
---|
4 | #include "evio_output.h" |
---|
5 | |
---|
6 | #include <fstream> |
---|
7 | |
---|
8 | void evio_output :: RecordAndClear(MOutputs* output, MBank mbank) |
---|
9 | { |
---|
10 | // Need to reorganize the vectors. |
---|
11 | // First index must be variable index |
---|
12 | // Second index must be hit number index |
---|
13 | |
---|
14 | vector<vector<double> > rawinfv; |
---|
15 | vector<vector<int> > dgtinfv; |
---|
16 | |
---|
17 | int NHITS = rawinfos.size(); |
---|
18 | if(NHITS) |
---|
19 | { |
---|
20 | int NVARSR = rawinfos[0].size(); |
---|
21 | int NVARSD = dgtinfos[0].size(); |
---|
22 | |
---|
23 | // cout << NHITS << " " << NVARSR << " " << NVARSD << endl; |
---|
24 | |
---|
25 | rawinfv.resize(NVARSR); |
---|
26 | dgtinfv.resize(NVARSD); |
---|
27 | |
---|
28 | for(int i=0; i<NVARSR; i++) rawinfv[i].resize(NHITS); |
---|
29 | for(int i=0; i<NVARSD; i++) dgtinfv[i].resize(NHITS); |
---|
30 | |
---|
31 | |
---|
32 | for(int nh=0; nh<NHITS; nh++) |
---|
33 | for(int iv=0; iv<NVARSR; iv++) |
---|
34 | rawinfv[iv][nh] = rawinfos[nh][iv]; |
---|
35 | |
---|
36 | for(int nh=0; nh<NHITS; nh++) |
---|
37 | for(int iv=0; iv<NVARSD; iv++) |
---|
38 | dgtinfv[iv][nh] = dgtinfos[nh][iv]; |
---|
39 | |
---|
40 | // inserting bitId in bankevent >> TAG=bankID NUM=1 << |
---|
41 | // hit size should be the same as bitId size |
---|
42 | if(NHITS != bitId.size()) cout << " Alarm! " << NHITS << " " << bitId.size() << endl; |
---|
43 | *bankevent << evioDOMNode::createEvioDOMNode(bankID, 1, bitId); |
---|
44 | |
---|
45 | // creating dgthitbank, inserting it in bankevent >> TAG=bankID NUM=100 << |
---|
46 | if(NVARSD) |
---|
47 | { |
---|
48 | evioDOMNodeP dgtbankhit = evioDOMNode::createEvioDOMNode(bankID, 100); |
---|
49 | *bankevent << dgtbankhit; |
---|
50 | |
---|
51 | for(int d=0; d<NVARSD; d++) |
---|
52 | *dgtbankhit << evioDOMNode::createEvioDOMNode(bankID, mbank.id[d+NVARSR], dgtinfv[d]); |
---|
53 | } |
---|
54 | |
---|
55 | // creating rawhitbank, inserting it in bankevent >> TAG=bankID NUM=200 << |
---|
56 | if(NVARSR) |
---|
57 | { |
---|
58 | evioDOMNodeP rawbankhit = evioDOMNode::createEvioDOMNode(bankID, 200); |
---|
59 | *bankevent << rawbankhit; |
---|
60 | |
---|
61 | for(int r=0; r<NVARSR; r++) |
---|
62 | *rawbankhit << evioDOMNode::createEvioDOMNode(bankID, mbank.id[r], rawinfv[r]); |
---|
63 | } |
---|
64 | } |
---|
65 | |
---|
66 | rawinfos.clear(); |
---|
67 | dgtinfos.clear(); |
---|
68 | bitId.clear(); |
---|
69 | } |
---|
70 | void evio_output :: WriteEvent(MOutputs* output) |
---|
71 | { |
---|
72 | output->pchan->write(*event); |
---|
73 | delete event; |
---|
74 | } |
---|
75 | |
---|
76 | void evio_output :: SetOutpHeader(int evtn, MOutputs* output) |
---|
77 | { |
---|
78 | event = new evioDOMTree(1, 0); |
---|
79 | |
---|
80 | // creating and inserting head bank >> TAG=1 NUM=1 << |
---|
81 | evioDOMNodeP head = evioDOMNode::createEvioDOMNode(1, 1); |
---|
82 | *event << evioDOMNode::createEvioDOMNode(1, 1, &evtn, 1); |
---|
83 | } |
---|
84 | |
---|
85 | |
---|
86 | |
---|
87 | void evio_output :: WriteGenerated(MOutputs* output, vector<MGeneratedParticle> MGP) |
---|
88 | { |
---|
89 | double MAXP = output->gemcOpt.args["NGENP"].arg; |
---|
90 | // creating and inserting generated particles bank >> TAG=10 NUM=10 << |
---|
91 | generatedp = evioDOMNode::createEvioDOMNode(10, 10); |
---|
92 | int NP = MGP.size(); |
---|
93 | vector< vector<double> > gensinfos; |
---|
94 | |
---|
95 | for(int i=0; i<MAXP && i<MGP.size(); i++) |
---|
96 | { |
---|
97 | vector<double> geninfo; |
---|
98 | geninfo.push_back((double) MGP[i].PID); |
---|
99 | geninfo.push_back(MGP[i].momentum.getX()/MeV); |
---|
100 | geninfo.push_back(MGP[i].momentum.getY()/MeV); |
---|
101 | geninfo.push_back(MGP[i].momentum.getZ()/MeV); |
---|
102 | geninfo.push_back(MGP[i].vertex.getX()/cm); |
---|
103 | geninfo.push_back(MGP[i].vertex.getY()/cm); |
---|
104 | geninfo.push_back(MGP[i].vertex.getZ()/cm); |
---|
105 | |
---|
106 | gensinfos.push_back(geninfo); |
---|
107 | } |
---|
108 | |
---|
109 | for(int p=0; p<gensinfos.size(); p++) *generatedp << evioDOMNode::createEvioDOMNode(1, 10*(p+1), gensinfos[p]); |
---|
110 | *event << generatedp; |
---|
111 | |
---|
112 | } |
---|
113 | |
---|
114 | |
---|
115 | |
---|
116 | void evio_output :: SetBankHeader(int bankid, string SDName, MOutputs* output) |
---|
117 | { |
---|
118 | bankID = bankid; |
---|
119 | bankevent = evioDOMNode::createEvioDOMNode(bankID, 0); |
---|
120 | *event << bankevent ; |
---|
121 | } |
---|
122 | |
---|
123 | void evio_output :: ProcessOutput(int bitcid, PH_output PHout, MOutputs* output, MBank mbank) |
---|
124 | { |
---|
125 | |
---|
126 | // check for double/int sizes consistency |
---|
127 | int nraws, ndigit; |
---|
128 | nraws=ndigit=0; |
---|
129 | for(int i=0; i<mbank.name.size(); i++) |
---|
130 | { |
---|
131 | if(mbank.type[i] == 0) ndigit++; |
---|
132 | if(mbank.type[i] == 1) nraws++; |
---|
133 | } |
---|
134 | if(PHout.raws.size() != nraws || PHout.dgtz.size() != ndigit) |
---|
135 | { |
---|
136 | cout << " Output does not match bank definition. This hit won't be written in the output stream." << endl; |
---|
137 | cout << " nraws size: " << nraws << " Output nraws: " << PHout.raws.size() << endl; |
---|
138 | cout << " ndgtz size: " << ndigit << " Output ndgt: " << PHout.dgtz.size() << endl; |
---|
139 | return; |
---|
140 | } |
---|
141 | |
---|
142 | vector<double> rawinf; |
---|
143 | vector<int> dgtinf; |
---|
144 | |
---|
145 | for(int r=0; r<nraws; r++) |
---|
146 | rawinf.push_back(PHout.raws[r]); |
---|
147 | |
---|
148 | for(int d=0; d<ndigit; d++) |
---|
149 | dgtinf.push_back(PHout.dgtz[d]); |
---|
150 | |
---|
151 | bitId.push_back(bitcid); |
---|
152 | rawinfos.push_back(rawinf); |
---|
153 | dgtinfos.push_back(dgtinf); |
---|
154 | } |
---|
155 | |
---|
156 | |
---|
157 | |
---|