Line data Source code
1 : #include "Plots.h"
2 : #include <iostream>
3 : #include <fstream>
4 : #include <math.h>
5 : using namespace std;
6 :
7 : namespace Tauolapp
8 : {
9 :
10 0 : Plots::Plots():
11 0 : m_incoming_pdg_id(1),
12 0 : m_cosTheta (-0.2),
13 0 : m_n_plot_points (1000)
14 0 : {
15 0 : }
16 :
17 : void Plots::SANCtest1(){
18 :
19 0 : cout<<"SANC plot 1 (short)..."<<endl;
20 :
21 : double smin = log(6.*6.)+0.0001;
22 : double smax = log(17000.*17000.);
23 0 : double step = (smax-smin)/(m_n_plot_points-1);
24 :
25 0 : ofstream f1,f2,f3;
26 0 : f1.open("f-sanc.txt");
27 0 : f2.open("f-born.txt");
28 0 : f3.open("f-plzap0.txt");
29 0 : for(int i=0; i<m_n_plot_points; i++)
30 : {
31 0 : double s = exp(smin+i*step);
32 :
33 : // Write SANC value
34 0 : t_pair.recalculateRij(m_incoming_pdg_id,15,s,m_cosTheta);
35 0 : f1<<sqrt(s)<<" "<<t_pair.m_R[0][3]<<endl;
36 :
37 : // Write Born-level value
38 : // (assuming table 11 is filled with born-level data)
39 0 : t_pair.recalculateRij(11,15,s,m_cosTheta);
40 0 : f2<<sqrt(s)<<" "<<t_pair.m_R[0][3]<<endl;
41 :
42 0 : int outgoing_pdg_id = 15;
43 :
44 : // Write PLZAP0 value
45 0 : double pz = 1-plzap0_(&m_incoming_pdg_id,&outgoing_pdg_id,&s, &m_cosTheta);
46 0 : t_pair.m_R[0][3]=2*pz-1;
47 0 : f3<<sqrt(s)<<" "<<t_pair.m_R[0][3]<<endl;
48 0 : }
49 0 : f1.close();
50 0 : f2.close();
51 0 : f3.close();
52 0 : }
53 :
54 : void Plots::SANCtest2(){
55 :
56 0 : cout<<"SANC plot 2 (short)..."<<endl;
57 :
58 : double smin = log(6.*6.)+0.0001;
59 : double smax = log(17000.*17000.);
60 0 : double step = (smax-smin)/(m_n_plot_points-1);
61 :
62 0 : ofstream f1,f2,f3;
63 0 : f1.open("f-w-single-point.txt");
64 0 : f2.open("f-w0-single-point.txt");
65 0 : f3.open("f-ww0-single-point.txt");
66 :
67 0 : for(int i=0; i<m_n_plot_points; i++){
68 :
69 0 : double s=exp(smin+i*step);
70 0 : t_pair.recalculateRij(m_incoming_pdg_id,15,s,m_cosTheta);
71 :
72 : // Write w, w0 and w/w0
73 0 : f1<<sqrt(s)<<" "<<Tauola::getEWwt()<<endl;
74 0 : f2<<sqrt(s)<<" "<<Tauola::getEWwt0()<<endl;
75 0 : f3<<sqrt(s)<<" "<<Tauola::getEWwt()/Tauola::getEWwt0()<<endl;
76 : }
77 0 : f1.close();
78 0 : f2.close();
79 0 : f3.close();
80 0 : }
81 :
82 : void Plots::SANCtest3(){
83 :
84 0 : cout<<"SANC plot 3 (long)..."<<endl;
85 :
86 : double smin = log(6.*6.)+0.0001;
87 : double smax = log(17000.*17000.);
88 0 : double step = (smax-smin)/(m_n_plot_points-1);
89 :
90 0 : ofstream f1;
91 0 : f1.open("f-err.txt");
92 : double costheta=-1.;
93 :
94 0 : for(int i=0; i<m_n_plot_points; i++){
95 :
96 : double buf=0.,err=0.;
97 :
98 0 : for(int j=0; j<m_n_plot_points; j++){
99 :
100 0 : double s = exp(smin+j*step);
101 :
102 : // Get value from SANC table
103 0 : t_pair.recalculateRij(m_incoming_pdg_id,15,s,costheta);
104 0 : buf = t_pair.m_R[0][3];
105 0 : t_pair.recalculateRij(11,15,s,costheta);
106 :
107 : // Calculate error between SANC and Born-level
108 0 : err += (buf-t_pair.m_R[0][3])*(buf-t_pair.m_R[0][3]);
109 : }
110 :
111 0 : f1<<costheta<<" "<<err/m_n_plot_points<<endl;
112 : err=0.;
113 0 : costheta+=2./m_n_plot_points;
114 : }
115 :
116 0 : f1.close();
117 0 : }
118 :
119 : void Plots::SANCtest4(){
120 :
121 0 : cout<<"SANC plot 4 (medium)..."<<endl;
122 :
123 : double smin = log(6.*6.);
124 : double smax = log(17000.*17000.);
125 0 : double step = (smax-smin)/(m_n_plot_points-1);
126 :
127 0 : ofstream f1,f2,f3;
128 0 : f1.open("f-cross.txt");
129 0 : f2.open("f-w.txt");
130 0 : f3.open("f-w0.txt");
131 :
132 0 : for(int i=0; i<m_n_plot_points; i++){
133 :
134 0 : double s = exp(smin+i*step);
135 : double sumEWwt = 0.;
136 : double sumEWwt0 = 0.;
137 : double costheta = -1.;
138 : int NCOS = 21;
139 :
140 : // Calculate sum of w/w0 for all cosTheta
141 0 : for(int j=0; j<NCOS; j++){
142 :
143 0 : costheta = -1. + 1.0/NCOS + j*2./NCOS;
144 :
145 0 : t_pair.recalculateRij(m_incoming_pdg_id,15,s,costheta);
146 :
147 0 : sumEWwt +=Tauola::getEWwt();
148 0 : sumEWwt0+=Tauola::getEWwt0();
149 : }
150 :
151 0 : f1<<sqrt(s)<<" "<<sumEWwt/sumEWwt0/m_n_plot_points<<endl;
152 0 : f2<<sqrt(s)<<" "<< 2./NCOS * sumEWwt <<endl;
153 0 : f3<<sqrt(s)<<" "<< 2./NCOS * sumEWwt0 <<endl;
154 : }
155 :
156 0 : f1.close();
157 0 : f2.close();
158 0 : f3.close();
159 0 : }
160 :
161 : void Plots::setSancVariables(int incoming, double cosTheta) {
162 0 : m_incoming_pdg_id = incoming;
163 0 : m_cosTheta = cosTheta;
164 0 : }
165 :
166 : } // namespace Tauolapp
|
