File indexing completed on 2023-10-25 09:34:31
0001 {
0002 gROOT->Reset();
0003 gROOT->SetStyle("Plain");
0004
0005 gStyle->SetOptStat(1111);
0006 gStyle->SetOptFit(111);
0007
0008 Float_t plsignal[50][73][5][5],plnoise[50][73][5][5];
0009 Float_t minsignal[50][73][5][5],minnoise[50][73][5][5];
0010 Float_t excluded_min[50][73][5][5];
0011
0012 FILE *Out1 = fopen("coefficients_8.9mln.txt", "w+");
0013
0014 for(Int_t k=1; k<50; k++)
0015 {
0016 for(Int_t i=1; i<73; i++)
0017 {
0018 for(Int_t j=1; j<5; j++)
0019 {
0020 for(Int_t l=1; l<5; l++)
0021 {
0022 excluded_min[k][i][j][l] = 0.;
0023 }
0024 }
0025 }
0026 }
0027
0028 std::string line;
0029
0030 std::ifstream in21( "HB_exclusion.txt" );
0031
0032 while( std::getline( in21, line)){
0033 istringstream linestream(line);
0034
0035 int eta,phi,dep;
0036 linestream>>dep>>eta>>phi;
0037 cout<<" Eta="<<eta<<endl;
0038 if( eta > 0 )
0039 {
0040 cout<<" eta > 0 "<<endl;
0041 } else
0042 {
0043 excluded_min[abs(eta)][phi][dep][1] = 1;
0044 }
0045 }
0046
0047
0048 cout<<" Read "<<endl;
0049
0050 std::ifstream in20( "var_noise_8.9mln.txt" );
0051
0052 while( std::getline( in20, line)){
0053 istringstream linestream(line);
0054 Float_t var,err;
0055 int subd,eta,phi,dep;
0056 linestream>>subd>>eta>>phi>>dep>>var;
0057
0058
0059 if( eta > 0 )
0060 {
0061 plnoise[eta][phi][dep][subd] = var;
0062 } else
0063 {
0064 minnoise[abs(eta)][phi][dep][subd] = var;
0065 }
0066 }
0067 cout<<" End of noise read "<<endl;
0068
0069 std::ifstream in21( "var_minbias_8.9mln.txt" );
0070
0071 while( std::getline( in21, line)){
0072 istringstream linestream(line);
0073 Float_t var,err;
0074 int subd,eta,phi,dep;
0075 linestream>>subd>>eta>>phi>>dep>>var>>err;
0076 if( eta > 0 )
0077 {
0078 plsignal[eta][phi][dep][subd] = var;
0079 } else
0080 {
0081 minsignal[abs(eta)][phi][dep][subd] = var;
0082 }
0083 }
0084
0085 cout<<" End of signal read "<<endl;
0086
0087
0088
0089
0090 Float_t plmean[50][5][5];
0091 Float_t minmean[50][5][5];
0092 Int_t numchan[50][5][5];
0093
0094 for(Int_t k=1; k<5; k++)
0095 {
0096 for(Int_t idep0=1; idep0<5; idep0++)
0097 {
0098 for(Int_t i=1; i<42; i++) {
0099
0100 plmean[i][idep0][k] = 0.;
0101 minmean[i][idep0][k] = 0.;
0102
0103 Int_t nch1 = 0;
0104 Int_t nchp = 0;
0105 Int_t nchm = 0;
0106
0107 for(Int_t j=1; j<73; j++){
0108 if(minsignal[i][j][idep0][k]>0.) nch1++;
0109 if(plsignal[i][j][idep0][k]>0.) nchp++;
0110 if(minsignal[i][j][idep0][k]>0.) nchm++;
0111 plmean[i][idep0][k] = plmean[i][idep0][k] + plsignal[i][j][idep0][k] - plnoise[i][j][idep0][k];
0112 minmean[i][idep0][k] = minmean[i][idep0][k] + minsignal[i][j][idep0][k] - minnoise[i][j][idep0][k];
0113 }
0114
0115 numchan[i][idep0][k] = nch1;
0116 if(nch1 == 0) continue;
0117 plmean[i][idep0][k] = plmean[i][idep0][k]/nchp;
0118 minmean[i][idep0][k] = minmean[i][idep0][k]/nchm;
0119 cout<<" k, idep0, i, nch1= "<<k<<" "<<idep0<<" "<<i<<" "<<nch1<<endl;
0120
0121 Float_t err = 0.00001;
0122
0123 if( k != 3 ) {
0124 for(Int_t j=1;j<73;j++) {
0125 if( plsignal[i][j][idep0][k] > 0.) {
0126 Float_t tt0;
0127 if( plsignal[i][j][idep0][k] - plnoise[i][j][idep0][k]>0.0000001 && plmean[i][idep0][k]>0.) {
0128 tt0 = sqrt(plmean[i][idep0][k]/(plsignal[i][j][idep0][k] - plnoise[i][j][idep0][k]));
0129 } else {
0130 tt0 = 1.;
0131 }
0132 fprintf(Out1,"%d %d %d %d %.5f %.5f\n",k,idep0,i,j,tt0,err);
0133 }
0134 }
0135 for(Int_t j=1;j<73;j++) {
0136 if( minsignal[i][j][idep0][k] > 0. ) {
0137 Float_t tt0;
0138 if(minsignal[i][j][idep0][k] - minnoise[i][j][idep0][k]>0.&& minmean[i][idep0][k]>0. ) {
0139 tt0 = sqrt(minmean[i][idep0][k]/(minsignal[i][j][idep0][k] - minnoise[i][j][idep0][k]));
0140 } else {
0141 tt0 = 1.;
0142 }
0143 Int_t ieta = -1*i;
0144 fprintf(Out1,"%d %d %d %d %.5f %.5f\n",k,idep0,ieta,j,tt0,err);
0145 }
0146 }
0147 }
0148 } else {
0149 for(Int_t j=1;j<73;j++){
0150 Float_t tt0=1.;
0151 Int_t ieta = i;
0152 fprintf(Out1,"%d %d %d %d %.5f %.5f\n",k,idep0,ieta,j,tt0,err);
0153 }
0154
0155 for(Int_t j=1;j<73;j++){
0156 Float_t tt0=1.;
0157 Int_t ieta = -1*i;
0158 fprintf(Out1,"%d %d %d %d %.5f %.5f\n",k,idep0,ieta,j,tt0,err);
0159 }
0160 }
0161 cout<<" End "<<endl;
0162 }
0163 }
0164 }
0165
0166 fclose(Out1);
0167 }