BOSS 7.1.2
BESIII Offline Software System
Loading...
Searching...
No Matches
TrackPool.cxx
Go to the documentation of this file.
4#include <string>
5
15 setBeamPosition(HepPoint3D(0.0,0.0,0.0));
16 setVBeamPosition(HepSymMatrix(3,0));
17 m_numberone = 0;
18 m_numbertwo = 0;
19}
20
21
22void TrackPool::AddTrack(const int number, const double mass,
23 const RecMdcTrack *trk) {
24 HepVector helix(5,0);
25 double error[15];
26 for(int i = 0; i < 5; i++)
27 helix[i] = trk->helix(i);
28 for(int i = 0; i < 15; i++)
29 error[i] = trk->err(i);
30 WTrackParameter wtrk(mass, helix, error);
31 setWTrackOrigin(wtrk);
32 setWTrackInfit(wtrk);
33 setWTrackList(number);
34 if(number != numberWTrack()-1) {
35 std::cout << "TrackPool: wrong track index" <<" "
36 <<number<<" , " <<numberWTrack()<< std::endl;
37 }
39 setMappositionA(m_numberone);
40 setMappositionB(m_numbertwo);
41 m_numberone = m_numberone + 4;
42}
43
44void TrackPool::AddTrack(const int number, const double mass,
45 const RecEmcShower *trk) {
46 //
47 //parameters: phi lambda mass energy
48 //
49 double ptrk = trk->energy();
50 double e = sqrt(ptrk*ptrk + mass * mass);
51 double the = trk->theta();
52 double phi = trk->phi();
53 HepLorentzVector p4(ptrk * sin(the) * cos(phi),
54 ptrk * sin(the) * sin(phi),
55 ptrk * cos(the),
56 e);
57 double dphi = trk->dphi();
58 double dthe = trk->dtheta();
59 double de = trk->dE();
60 double x = trk->x();
61 double y = trk->y();
62 double z = trk->z();
63 HepPoint3D x3 (x, y ,z);
64 WTrackParameter wtrk(x3, p4 ,dphi ,dthe, de);
65 HepSymMatrix Vpl = HepSymMatrix(2,0);
66 //=== get Vclus===
67 HepSymMatrix Vclus = HepSymMatrix (3,0);
68 Vclus = (wtrk.Ew()).sub(5,7);
69 double xpr = x - m_BeamPosition[0];
70 double ypr = y - m_BeamPosition[1];
71 double zpr = z - m_BeamPosition[2];
72 double Rpr = sqrt(xpr*xpr + ypr*ypr);
73 // === get jacobi ===
74 HepMatrix J(2,3,0);
75 J[0][0] = -ypr/(Rpr*Rpr);
76 J[0][1] = xpr/(Rpr*Rpr);
77 J[1][0] = -xpr * zpr/(Rpr*Rpr*Rpr);
78 J[1][1] = -ypr * zpr/(Rpr*Rpr*Rpr);
79 J[1][2] = 1/Rpr;
80 Vpl = Vclus.similarity(J) + m_VBeamPosition.similarity(J);
81 Vpl[0][1]=0;
82 // === get phipre, lambda===
83
84 double phipre = atan(ypr/xpr);
85
86 if(xpr<0){
87 phipre = atan(ypr/xpr) + 3.1415926;
88 }
89 double lambdapre = zpr/Rpr;
90
91
92 // === set p4 ===
93 double p0x = ptrk*cos(phipre)/sqrt(1 + lambdapre*lambdapre);
94 double p0y = ptrk*sin(phipre)/sqrt(1 + lambdapre*lambdapre);
95 double p0z = ptrk*lambdapre/sqrt(1 + lambdapre*lambdapre);
96 double p0e = e;
97
98
99 double p0ver = sqrt(p0x*p0x + p0y*p0y);
100
101
102 HepMatrix B(4,3,0);
103 B[0][0] = -p0y;
104 B[0][1] = -p0z * p0x * p0ver/(p0e * p0e);
105 B[0][2] = p0x/p0e;
106 B[1][0] = p0x;
107 B[1][1] = -p0z * p0y * p0ver/(p0e * p0e);
108 B[1][2] = p0y/p0e;
109 B[2][1] = p0ver * p0ver * p0ver/(p0e * p0e);
110 B[2][2] = p0z/p0e;
111 B[3][2] = 1;
112
113 HepSymMatrix Vple(3,0);
114 Vple[0][0] = Vpl[0][0];
115 Vple[1][1] = Vpl[1][1];
116 Vple[2][2] = de * de;
117
118 HepSymMatrix Vpxyze(4,0);
119 Vpxyze = Vple.similarity(B);
120
121 wtrk.setW(0,p0x);
122 wtrk.setW(1,p0y);
123 wtrk.setW(2,p0z);
124 wtrk.setW(3,p0e);
125
126 wtrk.setEw(Vpxyze);
127
128 HepSymMatrix Vplme(4,0);
129 Vplme[0][0] = Vpl[0][0];
130 Vplme[1][1] = Vpl[1][1];
131 Vplme[3][3] = de * de;
132 wtrk.setVplm(Vplme);
133
134 HepVector plmp(4 , 0);
135 plmp[0] = phipre;
136 plmp[1] = lambdapre;
137 plmp[2] = mass;
138 plmp[3] = e;
139 wtrk.setPlmp(plmp);
140
141
142 setWTrackOrigin(wtrk);
143 setWTrackInfit(wtrk);
144 setWTrackList(number);
145 if(number != numberWTrack()-1) {
146 std::cout << "TrackPool: wrong track index" <<" "
147 <<number<<" , " <<numberWTrack()<< std::endl;
148 }
149 GammaShape gtrk(p4,dphi,dthe,de);
150 setGammaShape(gtrk);
151 setGammaShapeList(number);
153 setMappositionA(m_numberone);
154 setMappositionB(m_numbertwo);
155 m_numberone = m_numberone + 4;
156}
157
158
159void TrackPool::AddMissTrack(const int number, const double mass,
160 const RecEmcShower *trk) {
161 //
162 //parameters: phi lambda mass ptrk
163 //
164 double ptrk = trk->energy();
165 double e = sqrt(ptrk*ptrk + mass * mass);
166 double the = trk->theta();
167 double phi = trk->phi();
168 HepLorentzVector p4( e* sin(the) * cos(phi),
169 e * sin(the) * sin(phi),
170 e * cos(the),
171 e);
172 double dphi = trk->dphi();
173 double dthe = trk->dtheta();
174 double de = 1E+6;
175 double x = trk->x();
176 double y = trk->y();
177 double z = trk->z();
178
179 HepPoint3D x3 (x, y ,z);
180 WTrackParameter wtrk(x3, p4 ,dphi ,dthe, de);
181 HepSymMatrix Vpe = HepSymMatrix(2,0);
182 //=== get Vclus===
183 HepSymMatrix Vclus = HepSymMatrix (3,0);
184 Vclus = (wtrk.Ew()).sub(5,7);
185 double xpr = x - m_BeamPosition[0];
186 double ypr = y - m_BeamPosition[1];
187 double zpr = z - m_BeamPosition[2];
188 double Rpr = sqrt(xpr*xpr + ypr*ypr);
189 // === get jacobi ===
190 HepMatrix J(2,3,0);
191 J[0][0] = -ypr/(Rpr*Rpr);
192 J[0][1] = xpr/(Rpr*Rpr);
193 J[1][0] = -xpr * zpr/(Rpr*Rpr*Rpr);
194 J[1][1] = -ypr * zpr/(Rpr*Rpr*Rpr);
195 J[1][2] = 1/Rpr;
196 Vpe = Vclus.similarity(J) + m_VBeamPosition.similarity(J);
197 Vpe[0][1]=0;
198
199 double phipre = atan(ypr/xpr);
200
201 if(xpr<0){
202 phipre = atan(ypr/xpr) + 3.1415926;
203 }
204 double lambdapre = zpr/Rpr;
205
206
207 HepVector plmp(4 , 0);
208 plmp[0] = phipre;
209 plmp[1] = lambdapre;
210 plmp[2] = mass;
211 plmp[3] = ptrk;
212 wtrk.setPlmp(plmp);
213
214 HepSymMatrix Vplm(3,0);
215 Vplm[0][0] = Vpe[0][0];
216 Vplm[1][1] = Vpe[1][1];
217 wtrk.setVplm(Vplm);
218
219
220
221 // === set p4 ===
222 double p0x = ptrk*cos(phipre)/sqrt(1 + lambdapre*lambdapre);
223 double p0y = ptrk*sin(phipre)/sqrt(1 + lambdapre*lambdapre);
224 double p0z = ptrk*lambdapre/sqrt(1 + lambdapre*lambdapre);
225 double p0e = e;
226
227 wtrk.setW(0,p0x);
228 wtrk.setW(1,p0y);
229 wtrk.setW(2,p0z);
230 wtrk.setW(3,p0e);
231
232 wtrk.setType(1);
233 setWTrackOrigin(wtrk);
234 setWTrackInfit(wtrk);
235 setWTrackList(number);
236 GammaShape gtrk(p4,dphi,dthe,de);
237 setGammaShape(gtrk);
238 setGammaShapeList(number);
240 setMappositionA(m_numberone);
241 setMappositionB(m_numbertwo);
242
243 m_numberone = m_numberone + 3;
244 m_numbertwo = m_numbertwo + 1;
245}
246
247
248
249
250void TrackPool::AddMissTrack(const int number, const RecEmcShower *trk) {
251 //
252 //parameters: phi lambda mass E
253 //
254
255 double mass = 0;
256 double ptrk = trk->energy();
257 double e = sqrt(ptrk*ptrk + mass * mass);
258 double the = trk->theta();
259 double phi = trk->phi();
260 HepLorentzVector p4( e* sin(the) * cos(phi),
261 e * sin(the) * sin(phi),
262 e * cos(the),
263 e);
264 double dphi = trk->dphi();
265 double dthe = trk->dtheta();
266 double de = 1E+6;
267 double x = trk->x();
268 double y = trk->y();
269 double z = trk->z();
270
271 HepPoint3D x3 (x, y ,z);
272 WTrackParameter wtrk(x3, p4 ,dphi ,dthe, de);
273 HepSymMatrix Vpe = HepSymMatrix(2,0);
274 //=== get Vclus===
275 HepSymMatrix Vclus = HepSymMatrix (3,0);
276 Vclus = (wtrk.Ew()).sub(5,7);
277 double xpr = x - m_BeamPosition[0];
278 double ypr = y - m_BeamPosition[1];
279 double zpr = z - m_BeamPosition[2];
280 double Rpr = sqrt(xpr*xpr + ypr*ypr);
281 // === get jacobi ===
282 HepMatrix J(2,3,0);
283 J[0][0] = -ypr/(Rpr*Rpr);
284 J[0][1] = xpr/(Rpr*Rpr);
285 J[1][0] = -xpr * zpr/(Rpr*Rpr*Rpr);
286 J[1][1] = -ypr * zpr/(Rpr*Rpr*Rpr);
287 J[1][2] = 1/Rpr;
288 Vpe = Vclus.similarity(J) + m_VBeamPosition.similarity(J);
289 Vpe[0][1]=0;
290 double phipre = atan(ypr/xpr);
291
292 if(xpr<0){
293 phipre = atan(ypr/xpr) + 3.1415926;
294 }
295 double lambdapre = zpr/Rpr;
296
297
298 HepVector plmp(4 , 0);
299 plmp[0] = phipre;
300 plmp[1] = lambdapre;
301 plmp[2] = mass;
302 plmp[3] = e;
303 wtrk.setPlmp(plmp);
304
305 HepSymMatrix Vplm(2,0);
306 Vplm[0][0] = Vpe[0][0];
307 Vplm[1][1] = Vpe[1][1];
308 wtrk.setVplm(Vplm);
309
310
311
312 // === set p4 ===
313 double p0x = ptrk*cos(phipre)/sqrt(1 + lambdapre*lambdapre);
314 double p0y = ptrk*sin(phipre)/sqrt(1 + lambdapre*lambdapre);
315 double p0z = ptrk*lambdapre/sqrt(1 + lambdapre*lambdapre);
316 double p0e = e;
317
318 wtrk.setW(0,p0x);
319 wtrk.setW(1,p0y);
320 wtrk.setW(2,p0z);
321 wtrk.setW(3,p0e);
322
323 wtrk.setType(1);
324 setWTrackOrigin(wtrk);
325 setWTrackInfit(wtrk);
326 setWTrackList(number);
327 GammaShape gtrk(p4,dphi,dthe,de);
328 setGammaShape(gtrk);
329 setGammaShapeList(number);
331 setMappositionA(m_numberone);
332 setMappositionB(m_numbertwo);
333
334 m_numberone = m_numberone + 2;
335 m_numbertwo = m_numbertwo + 2;
336}
337
338
339
340
341void TrackPool::AddMissTrack(const int number, const double mass, const HepLorentzVector p4) {
342 //
343 //parameters: mass px py pz
344 //
345
346 double dphi = 1E+6;
347 double dthe = 1E+6;
348 double dE = 1E+6;
349 WTrackParameter wtrk(p4, dphi, dthe, dE);
350 HepVector plmp(4, 0);
351 double phipre = atan(p4[1]/p4[0]);
352
353 if(p4[0]<0){
354 phipre = atan(p4[1]/p4[0]) + 3.1415926;
355 }
356 plmp[0] = phipre;
357 plmp[1] = wtrk.Lambda();
358 plmp[2] = mass;
359 plmp[3] = p4[3];
360 HepSymMatrix Vplm(3, 0);
361 wtrk.setPlmp(plmp);
362 wtrk.setVplm(Vplm);
363 wtrk.setType(1);
364 setWTrackOrigin(wtrk);
365 setWTrackInfit(wtrk);
366 setWTrackList(number);
368 setMappositionA(m_numberone);
369 setMappositionB(m_numbertwo);
370 m_numberone = m_numberone + 1;
371 m_numbertwo = m_numbertwo + 3;
372}
373
374
375void TrackPool::AddMissTrack(const int number, const double mass) {
376 //
377 //parameters: mass px py pz, but px,py,pz 's error matrix is set as 10e^6
378 //
379 WTrackParameter wtrk;
380 wtrk.setMass(mass);
381 HepVector w(7,0);
382 HepSymMatrix Ew(7,0);
383 w[0] = 0.2;
384 w[1] = 0.2;
385 w[2] = 0.2;
386 w[3] = sqrt(0.2*0.2*3 + mass*mass);
387 Ew[0][0] = 1E+6;
388 Ew[1][1] = 1E+6;
389 Ew[2][2] = 1E+6;
390 wtrk.setW(w);
391 wtrk.setEw(Ew);
392 setWTrackOrigin(wtrk);
393 setWTrackInfit(wtrk);
394 setWTrackList(number);
396 setMappositionA(m_numberone);
397 setMappositionB(m_numbertwo);
398 m_numberone = m_numberone + 4;
399
400}
401
402
403void TrackPool::AddMissTrack(const int number, HepLorentzVector p4) {
404 double dphi = 1E+3;
405 double dthe = 1E+3;
406 double dE = 1E+3;
407 WTrackParameter wtrk(p4, dphi, dthe, dE);
408 HepSymMatrix Ew = HepSymMatrix(7,0);
409 for (int i = 0; i < 7; i++) {
410 for (int j = 0; j < 7; j++) {
411 if(i==j) Ew[i][j] = 1E+6;
412 }
413 }
414 wtrk.setType(1);
415 wtrk.setEw(Ew);
416 setWTrackOrigin(wtrk);
417 setWTrackInfit(wtrk);
418 setWTrackList(number);
420 setMappositionA(m_numberone);
421 setMappositionB(m_numbertwo);
422 m_numbertwo = m_numbertwo + 4;
423}
424
425
426
427
428
429void TrackPool::AddTrack(const int number, WTrackParameter wtrk) {
430 setWTrackOrigin(wtrk);
431 setWTrackInfit(wtrk);
432 setWTrackList(number);
433 if(number != numberWTrack()-1) {
434 std::cout << "TrackPool: wrong track index" <<" "
435 <<number<<" , " <<numberWTrack()<< std::endl;
436 }
438 setMappositionA(m_numberone);
439 setMappositionB(m_numbertwo);
440 m_numberone = m_numberone + 4;
441}
442
443
444void TrackPool::AddTrackVertex(const int number, const double mass, const RecEmcShower *trk) {
445 double ptrk = trk->energy();
446 double e = sqrt(ptrk*ptrk + mass * mass);
447 double the = trk->theta();
448 double phi = trk->phi();
449 HepLorentzVector p4(ptrk * sin(the) * cos(phi),
450 ptrk * sin(the) * sin(phi),
451 ptrk * cos(the),
452 e);
453 double dphi = trk->dphi();
454 double dthe = trk->dtheta();
455 double de = trk->dE();
456 double x = trk->x();
457 double y = trk->y();
458 double z = trk->z();
459 HepPoint3D x3 (x, y ,z);
460 WTrackParameter wtrk(x3, p4 ,dphi ,dthe, de);
461 setWTrackOrigin(wtrk);
462 setWTrackInfit(wtrk);
463 setWTrackList(number);
464 if(number != numberWTrack()-1) {
465 std::cout << "TrackPool: wrong track index" <<" "
466 <<number<<" , " <<numberWTrack()<< std::endl;
467 }
468 GammaShape gtrk(p4,dphi,dthe,de);
469 setGammaShape(gtrk);
470 setGammaShapeList(number);
472 m_numbertwo = 0;
473 setMappositionA(m_numberone);
474 setMappositionB(m_numbertwo);
475
476 m_numberone = m_numberone + 4;
477 m_numbertwo = m_numbertwo + 3;
478}
479
480
481
482
483std::vector<int> TrackPool::AddList(int n1) {
484 std::vector<int> lis;
485 lis.clear();
486 lis.push_back(n1);
487 return lis;
488}
489
490std::vector<int> TrackPool::AddList(int n1, int n2) {
491 std::vector<int> lis;
492 lis.clear();
493 lis.push_back(n1);
494 lis.push_back(n2);
495 return lis;
496}
497
498std::vector<int> TrackPool::AddList(int n1, int n2, int n3) {
499 std::vector<int> lis;
500 lis.clear();
501 lis.push_back(n1);
502 lis.push_back(n2);
503 lis.push_back(n3);
504 return lis;
505}
506
507std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4) {
508 std::vector<int> lis;
509 lis.clear();
510 lis.push_back(n1);
511 lis.push_back(n2);
512 lis.push_back(n3);
513 lis.push_back(n4);
514 return lis;
515}
516
517std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5) {
518 std::vector<int> lis;
519 lis.clear();
520 lis.push_back(n1);
521 lis.push_back(n2);
522 lis.push_back(n3);
523 lis.push_back(n4);
524 lis.push_back(n5);
525 return lis;
526}
527
528std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6) {
529 std::vector<int> lis;
530 lis.clear();
531 lis.push_back(n1);
532 lis.push_back(n2);
533 lis.push_back(n3);
534 lis.push_back(n4);
535 lis.push_back(n5);
536 lis.push_back(n6);
537 return lis;
538}
539
540std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7) {
541 std::vector<int> lis;
542 lis.clear();
543 lis.push_back(n1);
544 lis.push_back(n2);
545 lis.push_back(n3);
546 lis.push_back(n4);
547 lis.push_back(n5);
548 lis.push_back(n6);
549 lis.push_back(n7);
550 return lis;
551}
552
553std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7, int n8) {
554 std::vector<int> lis;
555 lis.clear();
556 lis.push_back(n1);
557 lis.push_back(n2);
558 lis.push_back(n3);
559 lis.push_back(n4);
560 lis.push_back(n5);
561 lis.push_back(n6);
562 lis.push_back(n7);
563 lis.push_back(n8);
564 return lis;
565}
566
567std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7, int n8,
568 int n9) {
569 std::vector<int> lis;
570 lis.clear();
571 lis.push_back(n1);
572 lis.push_back(n2);
573 lis.push_back(n3);
574 lis.push_back(n4);
575 lis.push_back(n5);
576 lis.push_back(n6);
577 lis.push_back(n7);
578 lis.push_back(n8);
579 lis.push_back(n9);
580 return lis;
581}
582
583std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7, int n8,
584 int n9, int n10) {
585 std::vector<int> lis;
586 lis.clear();
587 lis.push_back(n1);
588 lis.push_back(n2);
589 lis.push_back(n3);
590 lis.push_back(n4);
591 lis.push_back(n5);
592 lis.push_back(n6);
593 lis.push_back(n7);
594 lis.push_back(n8);
595 lis.push_back(n9);
596 lis.push_back(n10);
597 return lis;
598}
599
600std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7, int n8,
601 int n9, int n10, int n11) {
602 std::vector<int> lis;
603 lis.clear();
604 lis.push_back(n1);
605 lis.push_back(n2);
606 lis.push_back(n3);
607 lis.push_back(n4);
608 lis.push_back(n5);
609 lis.push_back(n6);
610 lis.push_back(n7);
611 lis.push_back(n8);
612 lis.push_back(n9);
613 lis.push_back(n10);
614 lis.push_back(n11);
615 return lis;
616}
617
618std::vector<int> TrackPool::AddList(int n1, int n2, int n3, int n4, int n5, int n6, int n7, int n8,
619 int n9, int n10, int n11, int n12) {
620 std::vector<int> lis;
621 lis.clear();
622 lis.push_back(n1);
623 lis.push_back(n2);
624 lis.push_back(n3);
625 lis.push_back(n4);
626 lis.push_back(n5);
627 lis.push_back(n6);
628 lis.push_back(n7);
629 lis.push_back(n8);
630 lis.push_back(n9);
631 lis.push_back(n10);
632 lis.push_back(n11);
633 lis.push_back(n12);
634 return lis;
635}
636
637
double sin(const BesAngle a)
Definition BesAngle.h:210
double cos(const BesAngle a)
Definition BesAngle.h:213
double mass
Double_t x[10]
double w
HepGeom::Point3D< double > HepPoint3D
Definition Gam4pikp.cxx:37
int n2
Definition SD0Tag.cxx:55
int n1
Definition SD0Tag.cxx:54
double dphi() const
double theta() const
double phi() const
double x() const
double z() const
double energy() const
double dE() const
double dtheta() const
double y() const
const HepSymMatrix err() const
const HepVector helix() const
......
std::vector< int > AddList(int n1)
void setMappositionB(const int n)
Definition TrackPool.h:130
void setVBeamPosition(const HepSymMatrix VBeamPosition)
Definition TrackPool.h:144
int numberWTrack() const
Definition TrackPool.h:79
void AddTrack(const int number, const double mass, const RecMdcTrack *trk)
Definition TrackPool.cxx:22
void clearMapkinematic()
Definition TrackPool.h:124
void clearGammaShapeList()
Definition TrackPool.h:140
void clearWTrackOrigin()
Definition TrackPool.h:111
void clearMappositionB()
Definition TrackPool.h:126
void clearWTrackList()
Definition TrackPool.h:113
void AddTrackVertex(const int number, const double mass, const RecEmcShower *trk)
void clearMappositionA()
Definition TrackPool.h:125
void AddMissTrack(const int number, const double mass)
void setGammaShape(const int n, const GammaShape gammashape)
Definition TrackPool.h:135
void setWTrackList(const int n)
Definition TrackPool.h:109
void setGammaShapeList(const int n)
Definition TrackPool.h:137
void setMappositionA(const int n)
Definition TrackPool.h:129
void clearGammaShape()
Definition TrackPool.h:139
void clearWTrackInfit()
Definition TrackPool.h:112
void setBeamPosition(const HepPoint3D BeamPosition)
Definition TrackPool.h:143
void setMapkinematic(const int n)
Definition TrackPool.h:128
void setWTrackOrigin(const int n, const WTrackParameter wtrk)
Definition TrackPool.h:105
void setWTrackInfit(const int n, const WTrackParameter wtrk)
Definition TrackPool.h:106
void setEw(const HepSymMatrix &Ew)
void setMass(const double mass)
double Lambda() const
void setW(const HepVector &w)
void setType(const int type)
void setPlmp(const HepVector &plmp)
HepSymMatrix Ew() const
void setVplm(const HepSymMatrix &Vplm)
double y[1000]
float ptrk
double double double * p4
Definition qcdloop1.h:77