00001 #ifndef HOUGHVALIDUPDATE_H
00002 #define HOUGHVALIDUPDATE_H
00003 #include "GaudiKernel/Algorithm.h"
00004 #include "GaudiKernel/NTuple.h"
00005 #include "GaudiKernel/INTupleSvc.h"
00006 #include "MdcGeomSvc/IMdcGeomSvc.h"
00007 #include "MdcGeomSvc/MdcGeoWire.h"
00008 #include "MdcGeomSvc/MdcGeoLayer.h"
00009 #include "MdcPrintSvc/MdcPrintSvc.h"
00010 #include "RawDataProviderSvc/RawDataProviderSvc.h"
00011 #include <string>
00012 #include <vector>
00013 #include <iostream>
00014 #include <math.h>
00015 #include <map>
00016 #include <set>
00017 #include "TH2D.h"
00018
00019 #include "CLHEP/Alist/AList.h"
00020 #include "MdcGeom/MdcDetector.h"
00021 #include "BField/BField.h"
00022 #include "TrkBase/TrkRecoTrk.h"
00023 #include "TrkFitter/TrkContextEv.h"
00024 #include "MdcCalibFunSvc/MdcCalibFunSvc.h"
00025 #include "MdcUtilitySvc/MdcUtilitySvc.h"
00026 #include "MdcPrintSvc/MdcPrintSvc.h"
00027 #include "HepPDT/ParticleDataTable.hh"
00028 #include "MdcData/MdcHit.h"
00029 class HoughValidUpdate:public Algorithm {
00030 public:
00031 HoughValidUpdate(const std::string& name, ISvcLocator* pSvcLocator);
00032 StatusCode initialize();
00033 StatusCode execute();
00034 StatusCode finalize();
00035 StatusCode beginRun();
00036 typedef std::map<int, int> M2;
00037 typedef std::map<int, M2> M1;
00038 struct Mytrack{
00039 double pt;
00040 int charge;
00041 double phi;
00042 double x_cir;
00043 double y_cir;
00044 double r;
00045 bool use_in_tds;
00046 TrkRecoTrk* newTrack;
00047 vector<int> vec_hit_on_trk;
00048 };
00049 private:
00050 int GetMcInfo();
00051 int digiToHots(MdcDigiVec mdcDigiVec,TrkRecoTrk* newTrack,vector<bool> vec_truthHit,uint32_t getDigiFlag,vector<double> vec_flt_truth,map<int,double> hit_to_circle1,vector<MdcHit*>&);
00052 int digiToHots2(MdcDigiVec mdcDigiVec,TrkRecoTrk* newTrack,vector<bool> vec_truthHit,vector<double> vec_flt_truth,uint32_t getDigiFlag, vector<int> vec_slant,vector<double> vec_l_Calcu,vector<double> vec_l_axial,vector<double> vec_theta_ontrack,vector<int>& vec_hitbeforefit,map<int,double> hit_to_circle2,vector<MdcHit*>&,double );
00053 double CFtrans(double x,double y);
00054 int SlantId(int layer);
00055 int LeastSquare(int n,vector<double> n_x,vector<double> n_y, vector<int> n_slant,vector<int> n_layer,vector<bool> vec_truthHit,double &d0,double &phi0, double &omega,double& circleX, double& circleY,double& circleR);
00056 int Zposition(int n,vector<int> n_slant,double x_cirtemp,double y_cirtemp,double r_temp,vector<double> n_x_east,vector<double> n_x_west,vector<double> n_y_east,vector<double> n_y_west,vector<double> n_z_east,vector<double> n_z_west,vector<int> n_layer,vector<int> n_wire,vector<double> n_z,vector<double>& z_stereo_aver,vector<double>& l,vector<bool> vec_truthHit);
00057
00058 void Linefit(int z_stereoNum,vector<double> z_stereo_aver,vector<double> l, double& tanl,double& z0);
00059 void Multi_array(int binX,int binY);
00060
00061 void FillCells(TH2D *h1,int n,int binX,int binY,vector<double> vec_u,vector<double> vec_v,vector<int> vec_layer,vector<int> vec_wire,vector< vector <int> >& countij,vector< vector < vector<int> > >& vec_selectNum,vector< vector < vector<int> > >& vec_selectHit);
00062 void RhoTheta(int m_nHit,vector<double> vec_u,vector<double> vec_v,vector<double>& vec_rho,vector<double>& vec_theta,vector< vector<int> >& vec_hitNum, vector<int> vec_digitruth);
00063 void RhoThetaAxial(vector<int> vec_slant,int nHit,vector<double> vec_u,vector<double> vec_v,vector<double>& vec_rho,vector<double>& vec_theta,vector< vector<int> >& vec_hitNum,vector<int> vec_digitruth);
00064 void FillHist(TH2D *h1,vector<double> vec_u,vector<double> vec_v,int m_nHit,vector< vector < vector<int> > > &vec_selectNum);
00065 void FillRhoTheta(TH2D *h1 ,vector<double> vec_theta,vector<double> vec_rho);
00066 int zPosition(MdcHit& h, double t0, double x_cirtemp,double y_cirtemp, double r_temp,int digiId,double& delta_z,double& delta_l,double& l_Calcu_Left,double& l_Calcu_Right,double& z_Calcu_Left,double& z_Calcu_Right,double& z_Calcu,double& l_Calcu,int& i_q);
00067 void AmbigChoose(vector< vector< vector<double> > > point,int n_use,vector<int>& ambigList,double& tanl,double& z0);
00068 void Combine(TH2D *h1,int widthX,int widthY,vector< vector < vector<int> > > vec_selectNum,M1 &hit_combine,M2 &peak_combine_num,int &ntrack,vector<int> peakList1[]);
00069 void Combine_two_cell(TH2D *h1,int xPeakCell,int yPeakCell,int xAround,int yaround,vector< vector < vector<int> > > vec_selectNum,int m,M1 &hit_combine,M2 &peak_combine_num);
00070 int LeastLine(int ,double x[] ,double y[] ,double& ,double& ,double&);
00072
00073 int m_drawTuple;
00074 int m_useMcInfo;
00075 int m_par_use;
00076 int m_method;
00077 int m_debug;
00078 int m_data;
00079 int m_binx;
00080 int m_biny;
00081 int m_peakCellNumX;
00082 int m_peakCellNumY;
00083 int m_ndev1;
00084 int m_ndev2;
00085 double m_fpro;
00086 double m_hit_pro;
00087 std::string m_pdtFile;
00088 bool m_pickHits;
00089 double m_fltcut;
00090 double m_nHitcut_low;
00091 double m_nHitcut_high;
00092 int m_eventcut;
00093 double m_rhocut;
00094 int m_setpeak_fit;
00095 int m_mcpar;
00096 int m_fillTruth;
00097
00098
00099 double m_bunchT0;
00100 int t_eventNum;
00101 int t_runNum;
00102 double m_maxrho;
00103
00104 int binX;
00105 int binY;
00106
00107 double d0_mc;
00108 double phi0_mc;
00109 double omega_mc;
00110 double dz_mc;
00111 double tanl_mc;
00112 int track_fit;
00113
00114 uint32_t m_getDigiFlag;
00115 int m_maxMdcDigi;
00116 bool m_keepBadTdc;
00117 bool m_dropHot;
00118 bool m_keepUnmatch;
00119 int m_minMdcDigi;
00120 int m_q;
00121 double m_dropTrkDrCut;
00122 double m_dropTrkDzCut;
00123 double m_dropTrkPtCut;
00124 double m_dropTrkChi2Cut;
00125
00126 std::string m_configFile;
00127 const MdcDetector* m_gm;
00128 HepPDT::ParticleDataTable* m_particleTable;
00129 TrkContextEv* m_context;
00130 BField* m_bfield;
00131 const MdcCalibFunSvc* m_mdcCalibFunSvc;
00132
00133 RawDataProviderSvc* m_rawDataProviderSvc;
00134 MdcGeomSvc* m_mdcGeomSvc;
00135 MdcPrintSvc* m_mdcPrintSvc;
00136 int t_eventNo;
00137 IMagneticFieldSvc* m_pIMF;
00138 int m_pid;
00139 int nfailure;
00140 bool m_combineTracking;
00141 std::vector<float> m_helixHitsSigma;
00142 double m_qualityFactor;
00143 bool m_removeBadTrack;
00144 int m_dropTrkNhitCut;
00145 double m_dropTrkChi2NdfCut;
00146
00147
00148
00149 NTuple::Tuple* ntuplehit;
00150 NTuple::Item<int> m_eventNum;
00151 NTuple::Item<int> m_eventTime;
00152 NTuple::Item<int> m_runNum;
00153 NTuple::Item<int> m_nHit_Mc;
00154 NTuple::Item<int> m_nHitAxial_Mc;
00155 NTuple::Array<int> m_layer_Mc;
00156 NTuple::Array<int> m_cell_Mc;
00157
00158
00159 NTuple::Item<int> m_nTrkMC;
00160 NTuple::Array<double> m_pidTruth;
00161 NTuple::Array<double> m_costaTruth;
00162 NTuple::Array<double> m_ptTruth;
00163 NTuple::Array<double> m_pzTruth;
00164 NTuple::Array<double> m_pTruth;
00165 NTuple::Array<double> m_qTruth;
00166 NTuple::Array<double> m_drTruth;
00167 NTuple::Array<double> m_phi0Truth;
00168 NTuple::Array<double> m_omegaTruth;
00169 NTuple::Array<double> m_dzTruth;
00170 NTuple::Array<double> m_tanl_mc;
00171
00172 NTuple::Item<int> m_nHit;
00173 NTuple::Item<int> m_nHitDigi;
00174 NTuple::Item<int> m_nHitAxial;
00175 NTuple::Item<int> m_nHitStereo;
00176 NTuple::Array<int> m_layer;
00177 NTuple::Array<int> m_cell;
00178 NTuple::Array<double> m_x_east;
00179 NTuple::Array<double> m_y_east;
00180 NTuple::Array<double> m_x_west;
00181 NTuple::Array<double> m_y_west;
00182 NTuple::Array<double> m_z_east;
00183 NTuple::Array<double> m_z_west;
00184 NTuple::Array<double> m_x;
00185 NTuple::Array<double> m_y;
00186 NTuple::Array<double> m_z;
00187 NTuple::Array<double> m_u;
00188 NTuple::Array<double> m_v;
00189 NTuple::Array<double> m_p;
00190 NTuple::Array<int> m_slant;
00191 NTuple::Item<int> m_layer_max;
00192 NTuple::Array<int> m_layerNhit;
00193 NTuple::Array<double> m_ltruth;
00194 NTuple::Array<double> m_ztruth;
00195 NTuple::Array<double> m_postruth;
00196 NTuple::Array<int> m_digi_truth;
00197 NTuple::Item<int> m_e_delta;
00198
00199 NTuple::Item<int> m_nCross;
00200 NTuple::Array<double> m_rho;
00201 NTuple::Array<double> m_theta;
00202 NTuple::Item<int> m_xybinNum;
00203 NTuple::Array<int> m_xybin;
00204 NTuple::Item<double> m_xsigma;
00205 NTuple::Item<double> m_ysigma;
00206
00207
00208
00209
00210
00211 NTuple::Item<int> m_npeak_1;
00212 NTuple::Item<int> m_npeak_2;
00213 NTuple::Item<int> m_trackFailure;
00214 NTuple::Item<int> m_npeak_after_tracking;
00215 NTuple::Array<int> m_nHitSelect;
00216 NTuple::Array<int> m_nHitAxialSelect;
00217 NTuple::Array<int> m_nHitStereoSelect;
00218 NTuple::Array<int> m_axiallose;
00219 NTuple::Array<int> m_stereolose;
00220
00221 NTuple::Item<int> m_trackNum_Mc;
00222 NTuple::Item<int> m_trackNum;
00223 NTuple::Item<int> m_trackNum_fit;
00224 NTuple::Item<int> m_trackNum_tds;
00225 NTuple::Array<double> m_x_circle;
00226 NTuple::Array<double> m_y_circle;
00227 NTuple::Array<double> m_r;
00228 NTuple::Array<double> m_chis_pt;
00229 NTuple::Array<double> m_d0;
00230 NTuple::Array<double> m_phi0;
00231 NTuple::Array<double> m_omega;
00232 NTuple::Array<double> m_z0;
00233 NTuple::Array<double> m_tanl;
00234 NTuple::Array<double> m_pt;
00235 NTuple::Array<double> m_pt2;
00236 NTuple::Array<double> m_pz;
00237 NTuple::Array<double> m_pxyz;
00238
00239 NTuple::Array<int> m_nFitFailure;
00240 NTuple::Array<int> m_2d_nFit;
00241 NTuple::Array<int> m_2dFit;
00242 NTuple::Array<int> m_3d_nFit;
00243 NTuple::Array<int> m_3dFit;
00244 NTuple::Array<double> m_pt2_rec_final;
00245 NTuple::Array<double> m_p_rec_final;
00246 NTuple::Array<double> m_d0_final;
00247 NTuple::Array<double> m_phi0_final;
00248 NTuple::Array<double> m_omega_final;
00249 NTuple::Array<double> m_z0_final;
00250 NTuple::Array<double> m_tanl_final;
00251 NTuple::Array<double> m_chis_3d_final;
00252 NTuple::Array<int> m_Q;
00253
00254 NTuple::Array<double> m_pt2_rec_tds;
00255 NTuple::Array<double> m_p_rec_tds;
00256 NTuple::Array<double> m_d0_tds;
00257 NTuple::Array<double> m_phi0_tds;
00258 NTuple::Array<double> m_omega_tds;
00259 NTuple::Array<double> m_z0_tds;
00260 NTuple::Array<double> m_tanl_tds;
00261 NTuple::Array<double> m_chis_3d_tds;
00262 NTuple::Array<int> m_Q_tds;
00263
00264
00265 double m_truthPos[43][288][3];
00266
00267 NTuple::Array<double> m_dist_to_circle;
00268
00269
00270 NTuple::Item<int> m_nHitStereo_use;
00271 NTuple::Array<double> m_lCalcuLeft;
00272 NTuple::Array<double> m_lCalcuRight;
00273 NTuple::Array<double> m_zCalcuLeft;
00274 NTuple::Array<double> m_zCalcuRight;
00275 NTuple::Array<double> m_zCalcu;
00276 NTuple::Array<double> m_lCalcu;
00277 NTuple::Array<double> m_delta_z;
00278 NTuple::Array<double> m_delta_l;
00279 NTuple::Array<int> m_amb;
00280 NTuple::Item<int> m_ambig_no_match;
00281 NTuple::Item<double> m_tanl_Cal;
00282 NTuple::Item<double> m_z0_Cal;
00283 NTuple::Item<double> m_ambig_no_match_propotion;
00284
00285 NTuple::Item<int> m_nhitdis0;
00286 NTuple::Item<int> m_nhitdis1;
00287 NTuple::Item<int> m_nhitdis2;
00288 NTuple::Array<double> m_hit_dis0;
00289 NTuple::Array<double> m_hit_dis1;
00290 NTuple::Array<double> m_hit_dis2;
00291 NTuple::Array<double> m_hit_dis1_3d;
00292 NTuple::Array<double> m_hit_dis2_3d;
00293
00294
00295
00296
00297
00298
00299
00300
00301 NTuple::Item<bool> m_decay;
00302 NTuple::Item<int> m_outputTrk;
00303 NTuple::Array<int> m_hitOnTrk;
00304
00305 NTuple::Item<int> m_track_use_intrk1;
00306 NTuple::Item<int> m_noise_intrk1;
00307
00308
00309
00310 };
00311 bool compare(const HoughValidUpdate::Mytrack& a,const HoughValidUpdate::Mytrack& b);
00312 #endif
