Belle II Software  release-08-01-10
Sv123.cc
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 #include <ecl/digitization/algorithms.h>
9 #include <cmath>
10 using namespace std;
11 
12 
13 double Belle2::ECL::Sv123(double t, double t01, double tb1, double t02, double tb2, double td1, double ts1)
14 {
15 
16  double sv123 = 0.;
17  double dks0, dks1, dksm,
18  dw0, dw1, dwp, dwm, das1, dac1, das0, dac0, dzna, dksm2, ds, dd,
19  dcs0, dsn0, dzn0, td, ts, dr,
20  dcs0s, dsn0s, dcs0d, dsn0d, dcs1s, dsn1s, dcs1d, dsn1d;
21 
22 
23  if (t < 0.) return 0.;
24 
25  dr = (ts1 - td1) / td1;
26  if (abs(dr) >= 0.0000001) {
27  td = td1;
28  ts = ts1;
29  } else {
30  td = td1;
31  if (ts1 > td1) {
32  ts = td1 * 1.00001;
33  } else {
34  ts = td1 * 0.99999;
35  }
36  }
37 
38  dr = ((t01 - t02) * (t01 - t02) + (tb1 - tb2) * (tb1 - tb2)) / ((t01) * (t01) + (tb1) * (tb1));
39  dks0 = 1.0 / t01;
40  dks1 = 1.0 / t02;
41 
42  if (dr < 0.0000000001) {
43 
44  if (dks0 > dks1) {
45  dks0 = dks1 * 1.00001;
46  } else {
47  dks0 = dks1 * 0.99999;
48  }
49  }
50 
51 
52 
53 
54  dksm = dks1 - dks0;
55 
56  ds = 1. / ts;
57  dd = 1. / td;
58 
59  dw0 = 1. / tb1;
60  dw1 = 1. / tb2;
61  dwp = dw0 + dw1;
62  dwm = dw1 - dw0;
63 
64  dksm2 = dksm * dksm;
65 
66  dzna = (dksm2 + dwm * dwm) * (dksm2 + dwp * dwp);
67 
68 
69  das0 = dw1 * (dksm2 + dwp * dwm);
70  dac0 = -2 * dksm * dw0 * dw1;
71  das1 = dw0 * (dksm2 - dwp * dwm);
72  dac1 = -dac0;
73 
74 
75 
76 
77 
78  dsn0 = (ds - dks0);
79  dcs0 = -dw0;
80  dzn0 = dcs0 * dcs0 + dsn0 * dsn0;
81 
82  dsn0s = (dsn0 * das0 - dcs0 * dac0) / dzn0;
83  dcs0s = (dcs0 * das0 + dsn0 * dac0) / dzn0;
84 
85  dsn0 = (ds - dks1);
86  dcs0 = -dw1;
87  dzn0 = dcs0 * dcs0 + dsn0 * dsn0;
88 
89  dsn1s = (dsn0 * das1 - dcs0 * dac1) / dzn0;
90  dcs1s = (dcs0 * das1 + dsn0 * dac1) / dzn0;
91 
92 
93  dsn0 = (dd - dks0);
94  dcs0 = -dw0;
95  dzn0 = dcs0 * dcs0 + dsn0 * dsn0;
96 
97  dsn0d = (dsn0 * das0 - dcs0 * dac0) / dzn0;
98  dcs0d = (dcs0 * das0 + dsn0 * dac0) / dzn0;
99 
100  dsn0 = (dd - dks1);
101  dcs0 = -dw1;
102  dzn0 = dcs0 * dcs0 + dsn0 * dsn0;
103 
104  dsn1d = (dsn0 * das1 - dcs0 * dac1) / dzn0;
105  dcs1d = (dcs0 * das1 + dsn0 * dac1) / dzn0;
106 
107  //cppcheck dr = (ts - td) / td;
108 
109 
110 
111 
112  sv123 = ((((dsn0s - dsn0d) * sin(dw0 * t)
113  + (dcs0s - dcs0d) * cos(dw0 * t)) * exp(-t * dks0)
114  - (dcs0s + dcs1s) * exp(-t * ds) + (dcs0d + dcs1d) * exp(-t * dd)
115  + ((dsn1s - dsn1d) * sin(dw1 * t)
116  + (dcs1s - dcs1d) * cos(dw1 * t)) * exp(-t * dks1)) / dzna / (ts - td));
117 
118  return sv123;
119 
120 
121 }