File indexing completed on 2024-04-06 12:19:55
0001 #include <iostream>
0002 #include <cmath>
0003
0004 using std::cout;
0005 using std::endl;
0006
0007 struct etmiss_vec {
0008 unsigned mag;
0009 unsigned phi;
0010 };
0011 etmiss_vec calculate_etmiss_vec (const int ex, const int ey) ;
0012 unsigned phFine(unsigned ph) {
0013 unsigned r = ph % 18;
0014 return (r<9 ? r : 17-r);
0015 }
0016
0017 int main(int argc, char **argv)
0018 {
0019 const unsigned maxbin=9;
0020 unsigned histTru[maxbin]={0,0,0,0,0,0,0,0,0};
0021 unsigned histGct[maxbin]={0,0,0,0,0,0,0,0,0};
0022 for (int ex=-50; ex<50; ex++) { for (int ey=-50; ey<50; ey++) {
0023
0024 if ((ex*ex+ey*ey)<1000 && (ex*ex+ey*ey)>=500) {
0025 etmiss_vec met = calculate_etmiss_vec(ex,ey);
0026 unsigned phiF = phFine(met.phi);
0027 if (phiF<=maxbin) { histGct[phiF]++; } else { cout << " unexpected phi " << phiF << endl; }
0028 float phiVal = (atan2( (float) ey, (float) ex))/M_PI;
0029 if (phiVal<0) { phiVal += 2.0; }
0030 unsigned phiTrue = ((unsigned) (36.*phiVal));
0031 unsigned phiT = phFine(phiTrue);
0032 if (phiT<=maxbin) { histTru[phiT]++; } else { cout << " unexpected phi " << phiT << endl; }
0033 if (met.phi != phiTrue) {
0034 cout << "ex " << ex << " ey " << ey << " gct phi " << met.phi
0035 << " tru phi " << phiTrue << endl;
0036 }
0037 }
0038 } }
0039 cout << "Phi bins Gct"; for (int b=0; b<maxbin; b++) { cout << " " << histGct[b]; } cout << endl;
0040 cout << "Phi bins Tru"; for (int b=0; b<maxbin; b++) { cout << " " << histTru[b]; } cout << endl;
0041 return 0;
0042 }
0043
0044 etmiss_vec calculate_etmiss_vec (const int ex, const int ey)
0045 {
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059
0060 etmiss_vec result;
0061
0062 unsigned eneCoarse, phiCoarse;
0063 unsigned eneCorect, phiCorect;
0064
0065 const unsigned root2fact = 181;
0066 const unsigned corrFact[11] = {24, 39, 51, 60, 69, 77, 83, 89, 95, 101, 106};
0067 const unsigned corrDphi[11] = { 0, 1, 1, 2, 2, 3, 3, 3, 3, 4, 4};
0068
0069 bool s[3];
0070 unsigned Mx, My, Mw;
0071
0072 unsigned Dx, Dy;
0073 unsigned Fx, Fy;
0074 unsigned eFact;
0075
0076 unsigned b,phibin;
0077 bool midphi=false;
0078
0079
0080
0081 My = (ey<0 ? -ey : ey);
0082 Mx = (ex<0 ? -ex : ex);
0083 Mw = (((Mx+My)*root2fact)+0x80)>>8;
0084
0085 s[0] = (ey<0);
0086 s[1] = (ex<0);
0087 s[2] = (My>Mx);
0088
0089 phibin = 0; b = 0;
0090 for (int i=0; i<3; i++) {
0091 if (s[i]) { b=1-b;} phibin = 2*phibin + b;
0092 }
0093
0094 unsigned m=(My>Mx ? My : Mx);
0095 eneCoarse = (Mw>m ? Mw : m );
0096 phiCoarse = phibin*9;
0097
0098
0099
0100
0101
0102 for (eFact=0; eFact<10; eFact++) {
0103 Fx = Mx;
0104 Fy = My;
0105 Dx = (Mx*corrFact[eFact])>>8;
0106 Dy = (My*corrFact[eFact])>>8;
0107 if ((Dx>=Fy) || (Dy>=Fx)) {midphi=false; break;}
0108 if ((Fx+Dx)>=(Fy-Dy) && (Fy+Dy)>=(Fx-Dx)) {midphi=true; break;}
0109 }
0110 eneCorect = (eneCoarse*(128+eFact))>>7;
0111 if (midphi ^ (b==1)) {
0112 phiCorect = phiCoarse + 8 - corrDphi[eFact];
0113 } else {
0114 phiCorect = phiCoarse + corrDphi[eFact];
0115 }
0116
0117
0118
0119 result.mag = eneCorect;
0120 result.phi = phiCorect;
0121
0122 return result;
0123 }