00001 #include "StandardNtuple/NtpStRecord.h"
00002 #include "CandNtupleSR/NtpSRRecord.h"
00003 #include "CandNtupleSR/NtpSREvent.h"
00004 #include "CandNtupleSR/NtpSRTrack.h"
00005 #include "CandNtupleSR/NtpSRShower.h"
00006 #include "CandNtupleSR/NtpSRStrip.h"
00007 #include "NueAna/NueAnaTools/SntpHelpers.h"
00008 #include "AnalysisNtuples/ANtpDefaultValue.h"
00009 #include "VertexFinder/NtpVtxFinder/NtpVtxFinder.h"
00010 #include "HighHitVars.h"
00011 #include "HighHitVarsAna.h"
00012 #include "MessageService/MsgService.h"
00013 #include "TPad.h"
00014 #include "TLatex.h"
00015 #include "TNtuple.h"
00016 #include "TFile.h"
00017 #include "TMath.h"
00018 #include <cassert>
00019 #include <vector>
00020 #include <algorithm>
00021 #include <iostream>
00022 #include <fstream>
00023
00024 using namespace std;
00025
00026
00027 HighHitVarsAna::HighHitVarsAna(HighHitVars &hhv):
00028 fHighHitVars(hhv)
00029 {
00030 }
00031
00032
00033 HighHitVarsAna::~HighHitVarsAna()
00034 {}
00035
00036 void HighHitVarsAna::Analyze(int evtn, RecRecordImp<RecCandHeader> *srobj)
00037 {
00038 if(srobj==0){
00039 return;
00040 }
00041
00042 if(((dynamic_cast<NtpStRecord *>(srobj))==0)&&
00043 ((dynamic_cast<NtpSRRecord *>(srobj))==0)){
00044 return;
00045 }
00046
00047 NtpSRShower *shower = SntpHelpers::GetPrimaryShower(evtn,srobj);
00048 Analyze(srobj,shower);
00049
00050
00051 }
00052
00053 void HighHitVarsAna::Analyze(RecRecordImp<RecCandHeader> *srobj, NtpSRShower* ntpShower=0)
00054 {
00055
00056 fHighHitVars.Reset();
00057
00058 if(srobj==0){
00059 return;
00060 }
00061
00062 if (ntpShower==0) {
00063 return;
00064 }
00065
00066 if(((dynamic_cast<NtpStRecord *>(srobj))==0)&&
00067 ((dynamic_cast<NtpSRRecord *>(srobj))==0)){
00068 return;
00069 }
00070
00071
00072 float A1, A2, A3, A4, A5, A6;
00073 float A0, B0, C0, D0, E0, F0;
00074
00075 A0=0;
00076 B0=0;
00077 C0=0;
00078 D0=0;
00079 E0=0;
00080 F0=0;
00081
00082 int showerstrips;
00083 showerstrips=ntpShower->nstrip;
00084
00085 int j, k, l, m, n;
00086 j = k = l = m = n = 0;
00087
00088 for(int i=0;i<showerstrips;i++){
00089 Int_t index = ntpShower->stp[i];
00090 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00091 if(!strip){
00092 continue;
00093 }
00094
00096 A1 = evtstp0mip[index] + evtstp1mip[index];
00097
00098 if(A1>A0){
00099 fHighHitVars.high_hit_1_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00100 fHighHitVars.high_hit_1_mip = A1;
00101 fHighHitVars.high_hit_1_strip = strip->strip;
00102 fHighHitVars.high_hit_1_plane = strip->plane;
00103 fHighHitVars.high_hit_1_planeview = strip->planeview;
00104 fHighHitVars.high_hit_1_zpos = strip->z;
00105 fHighHitVars.high_hit_1_tpos = strip->tpos;
00106 j=i;
00107 }
00108
00109 if(A1<=A0){
00110 continue;
00111 }
00112
00113 A0 = fHighHitVars.high_hit_1_mip;
00114
00115 }
00116
00118
00119
00120 Float_t hit1;
00121 hit1=fHighHitVars.high_hit_1_mip;
00122
00123 for(int i=0;i<showerstrips;i++){
00124 Int_t index = ntpShower->stp[i];
00125 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00126 if(!strip){
00127 continue;
00128 }
00129
00130 if(showerstrips>1){
00131
00132 A2 = evtstp0mip[index] + evtstp1mip[index];
00133
00134 if(A2<=hit1&&A2>B0&&i!=j){
00135 fHighHitVars.high_hit_2_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00136 fHighHitVars.high_hit_2_mip = A2;
00137 fHighHitVars.high_hit_2_strip = strip->strip;
00138 fHighHitVars.high_hit_2_plane = strip->plane;
00139 fHighHitVars.high_hit_2_planeview = strip->planeview;
00140 fHighHitVars.high_hit_2_zpos = strip->z;
00141 fHighHitVars.high_hit_2_tpos = strip->tpos;
00142 k=i;
00143 }
00144 if(A2<B0){
00145 continue;
00146 }
00147
00148 B0 = fHighHitVars.high_hit_2_mip;
00149 }
00150
00151 if(showerstrips<=1){
00152 fHighHitVars.high_hit_2_ph = ANtpDefVal::kFloat;
00153 fHighHitVars.high_hit_2_mip = ANtpDefVal::kFloat;
00154 fHighHitVars.high_hit_2_strip = ANtpDefVal::kInt;
00155 fHighHitVars.high_hit_2_plane = ANtpDefVal::kInt;
00156 fHighHitVars.high_hit_2_planeview = ANtpDefVal::kInt;
00157 fHighHitVars.high_hit_2_zpos = ANtpDefVal::kFloat;
00158 fHighHitVars.high_hit_2_tpos = ANtpDefVal::kFloat;
00159 }
00160 }
00162
00163
00164 Float_t hit2;
00165 hit2=fHighHitVars.high_hit_2_mip;
00166
00167
00168 for(int i=0;i<showerstrips;i++){
00169 Int_t index = ntpShower->stp[i];
00170 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00171 if(!strip){
00172 continue;
00173 }
00174
00175
00176 if(showerstrips>2){
00177
00178 A3 = evtstp0mip[index] + evtstp1mip[index];
00179
00180 if(A3<=hit2&&A3>C0&&i!=j&&i!=k){
00181 fHighHitVars.high_hit_3_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00182 fHighHitVars.high_hit_3_mip = A3;
00183 fHighHitVars.high_hit_3_strip = strip->strip;
00184 fHighHitVars.high_hit_3_plane = strip->plane;
00185 fHighHitVars.high_hit_3_planeview = strip->planeview;
00186 fHighHitVars.high_hit_3_zpos = strip->z;
00187 fHighHitVars.high_hit_3_tpos = strip->tpos;
00188 l=i;
00189 }
00190 if(A3<C0){
00191 continue;
00192 }
00193
00194 C0 = fHighHitVars.high_hit_3_mip;
00195 }
00196
00197 if(showerstrips<=2){
00198 fHighHitVars.high_hit_3_ph = ANtpDefVal::kFloat;
00199 fHighHitVars.high_hit_3_mip = ANtpDefVal::kFloat;
00200 fHighHitVars.high_hit_3_strip = ANtpDefVal::kInt;
00201 fHighHitVars.high_hit_3_plane = ANtpDefVal::kInt;
00202 fHighHitVars.high_hit_3_planeview = ANtpDefVal::kInt;
00203 fHighHitVars.high_hit_3_zpos = ANtpDefVal::kFloat;
00204 fHighHitVars.high_hit_3_tpos = ANtpDefVal::kFloat;
00205 }
00206 }
00208
00209
00210 Float_t hit3;
00211 hit3=fHighHitVars.high_hit_3_mip;
00212
00213 for(int i=0;i<showerstrips;i++){
00214 Int_t index = ntpShower->stp[i];
00215 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00216 if(!strip){
00217 continue;
00218 }
00219
00220
00221 if(showerstrips>3){
00222
00223 A4 = evtstp0mip[index] + evtstp1mip[index];
00224
00225 if(A4<=hit3&&A4>D0&&i!=j&&i!=k&&i!=l){
00226 fHighHitVars.high_hit_4_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00227 fHighHitVars.high_hit_4_mip = A4;
00228 fHighHitVars.high_hit_4_strip = strip->strip;
00229 fHighHitVars.high_hit_4_plane = strip->plane;
00230 fHighHitVars.high_hit_4_planeview = strip->planeview;
00231 fHighHitVars.high_hit_4_zpos = strip->z;
00232 fHighHitVars.high_hit_4_tpos = strip->tpos;
00233 m=i;
00234 }
00235 if(A4<D0){
00236 continue;
00237 }
00238
00239 D0 = fHighHitVars.high_hit_4_mip;
00240 }
00241
00242 if(showerstrips<=3){
00243 fHighHitVars.high_hit_4_ph = ANtpDefVal::kFloat;
00244 fHighHitVars.high_hit_4_mip = ANtpDefVal::kFloat;
00245 fHighHitVars.high_hit_4_strip = ANtpDefVal::kInt;
00246 fHighHitVars.high_hit_4_plane = ANtpDefVal::kInt;
00247 fHighHitVars.high_hit_4_planeview = ANtpDefVal::kInt;
00248 fHighHitVars.high_hit_4_zpos = ANtpDefVal::kFloat;
00249 fHighHitVars.high_hit_4_tpos = ANtpDefVal::kFloat;
00250 }
00251 }
00253
00254 Float_t hit4;
00255 hit4=fHighHitVars.high_hit_4_mip;
00256
00257
00258 for(int i=0;i<showerstrips;i++){
00259 Int_t index = ntpShower->stp[i];
00260 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00261 if(!strip){
00262 continue;
00263 }
00264
00265 if(showerstrips>4){
00266
00267 A5 = evtstp0mip[index] + evtstp1mip[index];
00268
00269 if(A5<=hit4&&A5>E0&&i!=j&&i!=k&&i!=l&&i!=m){
00270 fHighHitVars.high_hit_5_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00271 fHighHitVars.high_hit_5_mip = A5;
00272 fHighHitVars.high_hit_5_strip = strip->strip;
00273 fHighHitVars.high_hit_5_plane = strip->plane;
00274 fHighHitVars.high_hit_5_planeview = strip->planeview;
00275 fHighHitVars.high_hit_5_zpos = strip->z;
00276 fHighHitVars.high_hit_5_tpos = strip->tpos;
00277 n=i;
00278 }
00279 if(A5<E0){
00280 continue;
00281 }
00282
00283 E0 = fHighHitVars.high_hit_5_mip;
00284 }
00285
00286 if(showerstrips<=4){
00287 fHighHitVars.high_hit_5_ph = ANtpDefVal::kFloat;
00288 fHighHitVars.high_hit_5_mip = ANtpDefVal::kFloat;
00289 fHighHitVars.high_hit_5_strip = ANtpDefVal::kInt;
00290 fHighHitVars.high_hit_5_plane = ANtpDefVal::kInt;
00291 fHighHitVars.high_hit_5_planeview = ANtpDefVal::kInt;
00292 fHighHitVars.high_hit_5_zpos = ANtpDefVal::kFloat;
00293 fHighHitVars.high_hit_5_tpos = ANtpDefVal::kFloat;
00294 }
00295 }
00297
00298 Float_t hit5;
00299 hit5=fHighHitVars.high_hit_5_mip;
00300
00301 for(int i=0;i<showerstrips;i++){
00302 Int_t index = ntpShower->stp[i];
00303 NtpSRStrip *strip = SntpHelpers::GetStrip(index,srobj);
00304 if(!strip){
00305 continue;
00306 }
00307
00308 if(showerstrips>5){
00309
00310 A6 = evtstp0mip[index] + evtstp1mip[index];
00311
00312 if(A6<=hit5&&A6>F0&&i!=j&&i!=k&&i!=l&&i!=m&&i!=n){
00313 fHighHitVars.high_hit_6_ph = strip->ph0.sigcor+strip->ph1.sigcor;
00314 fHighHitVars.high_hit_6_mip = A6;
00315 fHighHitVars.high_hit_6_strip = strip->strip;
00316 fHighHitVars.high_hit_6_plane = strip->plane;
00317 fHighHitVars.high_hit_6_planeview = strip->planeview;
00318 fHighHitVars.high_hit_6_zpos = strip->z;
00319 fHighHitVars.high_hit_6_tpos = strip->tpos;
00320 }
00321 if(A6<F0){
00322 continue;
00323 }
00324
00325 F0 = fHighHitVars.high_hit_6_mip;
00326 }
00327
00328 if(showerstrips<=5){
00329 fHighHitVars.high_hit_6_ph = ANtpDefVal::kFloat;
00330 fHighHitVars.high_hit_6_mip = ANtpDefVal::kFloat;
00331 fHighHitVars.high_hit_6_strip = ANtpDefVal::kInt;
00332 fHighHitVars.high_hit_6_plane = ANtpDefVal::kInt;
00333 fHighHitVars.high_hit_6_planeview = ANtpDefVal::kInt;
00334 fHighHitVars.high_hit_6_zpos = ANtpDefVal::kFloat;
00335 fHighHitVars.high_hit_6_tpos = ANtpDefVal::kFloat;
00336 }
00337 }
00339
00340
00341 fHighHitVars.showervtx_plane=ntpShower->vtx.plane;
00342 fHighHitVars.showervtx_u_pos=ntpShower->vtx.u;
00343 fHighHitVars.showervtx_v_pos=ntpShower->vtx.v;
00344 fHighHitVars.showerbeg_plane=ntpShower->plane.beg;
00345
00346
00347
00348
00349
00350 }
00351
