Reference Manual

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    EvtPlaneProducer
00004 // Class:      EvtPlaneProducer
00005 // 
00013 //
00014 // Original Author:  Sergey Petrushanko
00015 //         Created:  Fri Jul 11 10:05:00 2008
00016 // $Id: EvtPlaneProducer.cc,v 1.9 2010/11/03 19:47:15 yilmaz Exp $
00017 //
00018 //
00019 
00020 
00021 // system include files
00022 #include <memory>
00023 #include <iostream>
00024 #include <time.h>
00025 #include "TMath.h"
00026 
00027 // user include files
00028 #include "FWCore/Framework/interface/Frameworkfwd.h"
00029 #include "FWCore/Framework/interface/EDProducer.h"
00030 #include "FWCore/Framework/interface/Event.h"
00031 #include "FWCore/Framework/interface/EventSetup.h"
00032 
00033 #include "FWCore/Framework/interface/MakerMacros.h"
00034 #include "FWCore/Framework/interface/ESHandle.h"
00035 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00036 #include "FWCore/Utilities/interface/InputTag.h"
00037 
00038 #include "DataFormats/Candidate/interface/Candidate.h"
00039 #include "DataFormats/HeavyIonEvent/interface/EvtPlane.h"
00040 
00041 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
00042 #include "DataFormats/Common/interface/EDProduct.h"
00043 #include "DataFormats/Common/interface/Ref.h"
00044 
00045 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
00046 #include "HepMC/GenEvent.h"
00047 #include "HepMC/GenParticle.h"
00048 
00049 #include "DataFormats/Common/interface/Handle.h"
00050 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
00051 
00052 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00053 #include "DataFormats/TrackReco/interface/Track.h"
00054 
00055 #include <cstdlib>
00056 using std::rand;
00057 
00058 using namespace std;
00059 
00060 //
00061 // class decleration
00062 //
00063 
00064 class EvtPlaneProducer : public edm::EDProducer {
00065    public:
00066       explicit EvtPlaneProducer(const edm::ParameterSet&);
00067       ~EvtPlaneProducer();
00068 
00069    private:
00070       virtual void beginJob() ;
00071       virtual void produce(edm::Event&, const edm::EventSetup&);
00072       virtual void endJob() ;
00073       
00074       // ----------member data ---------------------------
00075 
00076   bool useECAL_;
00077   bool useHCAL_;
00078   bool useTrackMidEta_;
00079   bool useTrackPosEta_;
00080   bool useTrackNegEta_;
00081   bool useTrackEta_;
00082   bool useTrackOdd_;
00083   bool useTrackEven_;
00084   bool useTrackPosCh_;
00085   bool useTrackNegCh_;
00086   bool useTrackPosEtaGap_;
00087   bool useTrackNegEtaGap_;
00088 
00089 
00090   edm::InputTag mcSrc_;
00091   edm::InputTag trackSrc_;
00092   edm::InputTag towerSrc_;
00093 
00094 
00095 };
00096 
00097 //
00098 // constants, enums and typedefs
00099 //
00100 
00101 
00102 //
00103 // static data member definitions
00104 //
00105 
00106 //
00107 // constructors and destructor
00108 //
00109 EvtPlaneProducer::EvtPlaneProducer(const edm::ParameterSet& iConfig)
00110 {
00111 
00112 int k=100;
00113 int iseed=k;
00114 
00115   srand(iseed);
00116 
00117    //register your products
00118   
00119   useECAL_ = iConfig.getUntrackedParameter<bool>("useECAL",true);
00120   useHCAL_ = iConfig.getUntrackedParameter<bool>("useHCAL",true);
00121   useTrackMidEta_ = iConfig.getUntrackedParameter<bool>("useTrackMidEta",true);
00122   useTrackPosEta_ = iConfig.getUntrackedParameter<bool>("useTrackPosEta",true);
00123   useTrackNegEta_ = iConfig.getUntrackedParameter<bool>("useTrackNegEta",true);
00124   useTrackEta_ = iConfig.getUntrackedParameter<bool>("useTrackEta",true);
00125   useTrackOdd_ = iConfig.getUntrackedParameter<bool>("useTrackOdd",true);
00126   useTrackEven_ = iConfig.getUntrackedParameter<bool>("useTrackEven",true);
00127   useTrackPosCh_ = iConfig.getUntrackedParameter<bool>("useTrackPosCh",true);
00128   useTrackNegCh_ = iConfig.getUntrackedParameter<bool>("useTrackNegCh",true);
00129   useTrackPosEtaGap_ = iConfig.getUntrackedParameter<bool>("useTrackPosEtaGap",true);
00130   useTrackNegEtaGap_ = iConfig.getUntrackedParameter<bool>("useTrackNegEtaGap",true);
00131 
00132 
00133   mcSrc_ = iConfig.getUntrackedParameter<edm::InputTag>("mc",edm::InputTag("generator"));
00134   trackSrc_ = iConfig.getUntrackedParameter<edm::InputTag>("tracks",edm::InputTag("hiSelectedTracks"));
00135   towerSrc_ = iConfig.getUntrackedParameter<edm::InputTag>("towers",edm::InputTag("towerMaker"));
00136 
00137   produces<reco::EvtPlaneCollection>("recoLevel");
00138 
00139 }
00140 
00141 
00142 EvtPlaneProducer::~EvtPlaneProducer()
00143 {
00144   
00145   // do anything here that needs to be done at desctruction time
00146   // (e.g. close files, deallocate resources etc.)
00147 
00148 }
00149 
00150 
00151 //
00152 // member functions
00153 //
00154 
00155 // ------------ method called to produce the data  ------------
00156 void
00157 EvtPlaneProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
00158 {
00159   using namespace edm;
00160   using namespace reco;
00161   using namespace HepMC;
00162 
00163   //calorimetry part
00164 
00165       double ugol[24], ugol2[24];
00166       double etugol[24], etugol2[24];
00167       double tower_eta, tower_phi;
00168       double tower_energy, tower_energy_e, tower_energy_h;
00169       double tower_energyet, tower_energyet_e, tower_energyet_h;
00170       double s1t, s2t, s1e, s2e, s1h, s2h, s1t1, s2t1;
00171       double ets1t, ets2t, ets1e, ets2e, ets1h, ets2h, ets1t1, ets2t1;
00172 //      double pi = 3.14159;
00173       double s1[24], s2[24];
00174       double ets1[24], ets2[24];
00175       
00176 //      double TEnergy[144], TPhi[144];
00177 //      int numb;       
00178 
00179 //      double planeA     =  0.;
00180 
00181 //       cout << endl << "  Start of the event plane determination." << endl;
00182 
00183        for(int j=0;j<24;j++) {
00184         s1[j]  = 0.;
00185         s2[j]  = 0.;
00186         ets1[j]  = 0.;
00187         ets2[j]  = 0.;
00188        }
00189       
00190 //       for(int l=0;l<144;l++) {
00191 //        TEnergy[l]  = 0.;
00192 //        TPhi[l]  = 0.;
00193 //       }
00194 
00195       Handle<CaloTowerCollection> calotower;
00196       iEvent.getByLabel(towerSrc_,calotower);
00197       
00198       if(!calotower.isValid()){
00199         cout << "Error! Can't get calotower product!" << endl;
00200        return ;
00201       }
00202 
00203         for (CaloTowerCollection::const_iterator j = calotower->begin();j !=calotower->end(); j++) {
00204 
00205 //        cout << *j << std::endl;
00206 //        cout << "ENERGY HAD " << j->hadEnergy()<< " ENERGY EM " <<j->emEnergy() 
00207 //        << " ETA " <<j->eta() << " PHI " <<j->phi() << std::endl;
00208     
00209         tower_eta        = j->eta();
00210         tower_phi        = j->phi();
00211         tower_energy_e   = j->emEnergy();
00212         tower_energy_h   = j->hadEnergy();
00213         tower_energy     = tower_energy_e + tower_energy_h;
00214         tower_energyet_e   = j->emEt();
00215         tower_energyet_h   = j->hadEt();
00216         tower_energyet     = tower_energyet_e + tower_energyet_h;
00217         
00218         s1t = tower_energy*sin(2.*tower_phi);
00219         s2t = tower_energy*cos(2.*tower_phi);
00220         s1e = tower_energy_e*sin(2.*tower_phi);
00221         s2e = tower_energy_e*cos(2.*tower_phi);
00222         s1h = tower_energy_h*sin(2.*tower_phi);
00223         s2h = tower_energy_h*cos(2.*tower_phi);
00224 
00225         s1t1 = tower_energy*sin(tower_phi);
00226         s2t1 = tower_energy*cos(tower_phi);
00227 
00228         ets1t = tower_energyet*sin(2.*tower_phi);
00229         ets2t = tower_energyet*cos(2.*tower_phi);
00230         ets1e = tower_energyet_e*sin(2.*tower_phi);
00231         ets2e = tower_energyet_e*cos(2.*tower_phi);
00232         ets1h = tower_energyet_h*sin(2.*tower_phi);
00233         ets2h = tower_energyet_h*cos(2.*tower_phi);
00234 
00235         ets1t1 = tower_energyet*sin(tower_phi);
00236         ets2t1 = tower_energyet*cos(tower_phi);
00237 
00238 // ENERGY
00239           
00240            s1[12] += s1t;
00241            s2[12] += s2t;
00242           
00243            ets1[12] += ets1t;
00244            ets2[12] += ets2t;
00245 
00246           if (tower_eta>0.25)
00247           {
00248            s1[19] +=  s1t;
00249            s2[19] +=  s2t;        
00250           }
00251           
00252           if (tower_eta<-0.25)
00253           {
00254            s1[20] +=  s1t;
00255            s2[20] +=  s2t;        
00256           }     
00257 
00258          if (fabs(tower_eta)>3.){ 
00259           if (tower_eta>3.) {
00260            s1[9] += s1t;
00261            s2[9] += s2t;
00262            s1[21] += s1t1;
00263            s2[21] += s2t1;
00264            
00265           }
00266           if (tower_eta<-3.) {
00267            s1[10] += s1t;
00268            s2[10] += s2t;
00269            s1[22] += s1t1;
00270            s2[22] += s2t1;
00271           }
00272            s1[11] += s1t;
00273            s2[11] += s2t;
00274            s1[23] += s1t1;
00275            s2[23] += s2t1;
00276          }
00277 
00278          if (fabs(tower_eta)<3.){
00279 
00280           s1[0] +=  s1t;
00281           s2[0] +=  s2t;
00282           s1[3] +=  s1h;
00283           s2[3] +=  s2h;
00284           s1[6] +=  s1e;
00285           s2[6] +=  s2e;
00286           
00287           if (tower_eta>0.25)
00288           {
00289            s1[13] +=  s1e;
00290            s2[13] +=  s2e;
00291            s1[15] +=  s1h;
00292            s2[15] +=  s2h;
00293            s1[17] +=  s1t;
00294            s2[17] +=  s2t;        
00295           }
00296           
00297           if (tower_eta<-0.25)
00298           {
00299            s1[14] +=  s1e;
00300            s2[14] +=  s2e;
00301            s1[16] +=  s1h;
00302            s2[16] +=  s2h;
00303            s1[18] +=  s1t;
00304            s2[18] +=  s2t;        
00305           }
00306                   
00307          if (fabs(tower_eta)>1.5) {
00308           s1[2] +=  s1t;
00309           s2[2] +=  s2t;
00310           s1[5] +=  s1h;
00311           s2[5] +=  s2h;
00312           s1[8] +=  s1e;
00313           s2[8] +=  s2e;
00314          }
00315          }
00316          
00317          if (fabs(tower_eta)<1.5){
00318           s1[1] +=  s1t;
00319           s2[1] +=  s2t;
00320           s1[4] +=  s1h;
00321           s2[4] +=  s2h;
00322           s1[7] +=  s1e;
00323           s2[7] +=  s2e;
00324          }
00325          
00326           if (tower_eta>0.25)
00327           {
00328            ets1[19] +=  ets1t;
00329            ets2[19] +=  ets2t;    
00330           }
00331           
00332           if (tower_eta<-0.25)
00333           {
00334            ets1[20] +=  ets1t;
00335            ets2[20] +=  ets2t;    
00336           }     
00337 
00338          if (fabs(tower_eta)>3.){ 
00339           if (tower_eta>3.) {
00340            ets1[9] += ets1t;
00341            ets2[9] += ets2t;
00342            ets1[21] += ets1t1;
00343            ets2[21] += ets2t1;
00344           }
00345           if (tower_eta<-3.) {
00346            ets1[10] += ets1t;
00347            ets2[10] += ets2t;
00348            ets1[22] += ets1t1;
00349            ets2[22] += ets2t1;
00350           }
00351            ets1[11] += ets1t;
00352            ets2[11] += ets2t;
00353            ets1[23] += ets1t1;
00354            ets2[23] += ets2t1;
00355          }
00356 
00357          if (fabs(tower_eta)<3.){
00358 
00359           ets1[0] +=  ets1t;
00360           ets2[0] +=  ets2t;
00361           ets1[3] +=  ets1h;
00362           ets2[3] +=  ets2h;
00363           ets1[6] +=  ets1e;
00364           ets2[6] +=  ets2e;
00365           
00366           if (tower_eta>0.25)
00367           {
00368            ets1[13] +=  ets1e;
00369            ets2[13] +=  ets2e;
00370            ets1[15] +=  ets1h;
00371            ets2[15] +=  ets2h;
00372            ets1[17] +=  ets1t;
00373            ets2[17] +=  ets2t;    
00374           }
00375           
00376           if (tower_eta<-0.25)
00377           {
00378            ets1[14] +=  ets1e;
00379            ets2[14] +=  ets2e;
00380            ets1[16] +=  ets1h;
00381            ets2[16] +=  ets2h;
00382            ets1[18] +=  ets1t;
00383            ets2[18] +=  ets2t;    
00384           }
00385                   
00386          if (fabs(tower_eta)>1.5) {
00387           ets1[2] +=  ets1t;
00388           ets2[2] +=  ets2t;
00389           ets1[5] +=  ets1h;
00390           ets2[5] +=  ets2h;
00391           ets1[8] +=  ets1e;
00392           ets2[8] +=  ets2e;
00393          }
00394          }
00395          
00396          if (fabs(tower_eta)<1.5){
00397           ets1[1] +=  ets1t;
00398           ets2[1] +=  ets2t;
00399           ets1[4] +=  ets1h;
00400           ets2[4] +=  ets2h;
00401           ets1[7] +=  ets1e;
00402           ets2[7] +=  ets2e;
00403          }       
00404          
00405         }
00406         
00407       for(int j1=0;j1<24;j1++) {
00408  
00409        if (s2[j1]==0.) {ugol[j1]=0.;}
00410        else {ugol[j1] = 0.5*atan2(s1[j1],s2[j1]);}
00411        
00412        if ((j1==21) || (j1==22) || (j1==23))
00413       {
00414        if (s2[j1]==0.) {ugol[j1]=0.;}
00415        else {ugol[j1] = atan2(s1[j1],s2[j1]);}       
00416        }
00417        
00418        ugol2[j1] = ugol[j1];
00419        
00420       }
00421 
00422 /*
00423        cout <<  endl << "   Azimuthal angle of reaction plane (with maximum), energy" << endl
00424        << "HCAL+ECAL (b+e)   " << ugol2[0] << endl
00425        << "HCAL+ECAL (b)     " << ugol2[1] << endl
00426        << "HCAL+ECAL (e)     " << ugol2[2] << endl
00427        << "HCAL      (b+e)   " << ugol2[3] << endl
00428        << "HCAL      (b)     " << ugol2[4] << endl
00429        << "HCAL      (e)     " << ugol2[5] << endl
00430        << "ECAL      (b+e)   " << ugol2[6] << endl
00431        << "ECAL      (b)     " << ugol2[7] << endl
00432        << "ECAL      (e)     " << ugol2[8] << endl  
00433        << "HF       (plus)   " << ugol2[9] << endl
00434        << "HF       (minus)  " << ugol2[10] << endl
00435        << "HF       (p+m)    " << ugol2[11] << endl
00436        << "HCAL+ECAL+HF      " << ugol2[12] << endl
00437        << "ECAL      (plus)  " << ugol2[13] << endl
00438        << "ECAL      (minus) " << ugol2[14] << endl
00439        << "HCAL      (plus)  " << ugol2[15] << endl
00440        << "HCAL      (minus) " << ugol2[16] << endl
00441        << "HCAL+ECAL (plus)  " << ugol2[17] << endl
00442        << "HCAL+ECAL (minus) " << ugol2[18] << endl
00443        << "HCAL+ECAL+HF (p)  " << ugol2[19] << endl
00444        << "HCAL+ECAL+HF (m)  " << ugol2[20] << endl
00445        << endl;
00446 */
00447 
00448 
00449 // ET ENERGY
00450 
00451         
00452       for(int j1=0;j1<24;j1++) {
00453  
00454        if (ets2[j1]==0.) {etugol[j1]=0.;}
00455        else {etugol[j1] = 0.5*atan2(ets1[j1],ets2[j1]);}
00456        
00457        if ((j1==21) || (j1==22) || (j1==23))
00458       {
00459        if (ets2[j1]==0.) {etugol[j1]=0.;}
00460        else {etugol[j1] = atan2(ets1[j1],ets2[j1]);}       
00461        }
00462               
00463        etugol2[j1] = etugol[j1];
00464        
00465       }
00466 
00467 /*
00468        cout <<  endl << "   Azimuthal angle of reaction plane (with maximum), t energy" << endl
00469        << "HCAL+ECAL (b+e)   " << etugol2[0] << endl
00470        << "HCAL+ECAL (b)     " << etugol2[1] << endl
00471        << "HCAL+ECAL (e)     " << etugol2[2] << endl
00472        << "HCAL      (b+e)   " << etugol2[3] << endl
00473        << "HCAL      (b)     " << etugol2[4] << endl
00474        << "HCAL      (e)     " << etugol2[5] << endl
00475        << "ECAL      (b+e)   " << etugol2[6] << endl
00476        << "ECAL      (b)     " << etugol2[7] << endl
00477        << "ECAL      (e)     " << etugol2[8] << endl  
00478        << "HF       (plus)   " << etugol2[9] << endl
00479        << "HF       (minus)  " << etugol2[10] << endl
00480        << "HF       (p+m)    " << etugol2[11] << endl
00481        << "HCAL+ECAL+HF      " << etugol2[12] << endl
00482        << "ECAL      (plus)  " << etugol2[13] << endl
00483        << "ECAL      (minus) " << etugol2[14] << endl
00484        << "HCAL      (plus)  " << etugol2[15] << endl
00485        << "HCAL      (minus) " << etugol2[16] << endl
00486        << "HCAL+ECAL (plus)  " << etugol2[17] << endl
00487        << "HCAL+ECAL (minus) " << etugol2[18] << endl
00488        << "HCAL+ECAL+HF (p)  " << etugol2[19] << endl
00489        << "HCAL+ECAL+HF (m)  " << etugol2[20] << endl
00490        << endl;
00491 
00492 
00493    const GenEvent *evt;
00494 
00495    Handle<HepMCProduct> mc;
00496    iEvent.getByLabel(mcSrc_,mc);
00497    evt = mc->GetEvent();
00498 
00499    const HeavyIon* hi = evt->heavy_ion();
00500    double impactP = hi->impact_parameter();
00501    double phi0 = hi->event_plane_angle();
00502 
00503    cout<<"  Generator Ugol " << phi0 <<  " Impact parameter  " << impactP << endl;
00504    
00505 */
00506 
00507 
00508 //Tracking part
00509      
00510       double track_eta;
00511       double track_phi;
00512       double track_pt;
00513       double track_charge;
00514 
00515    double trackPsi_eta_mid;
00516    double trackPsi_eta_pos;
00517    double trackPsi_eta_neg;
00518    double trackPsi_eta_posgap;
00519    double trackPsi_eta_neggap;
00520    double trackPsi_odd;
00521    double trackPsi_even;
00522    double trackPsi_eta;
00523    double trackPsi_PosCh;
00524    double trackPsi_NegCh;
00525 
00526 
00527    Handle<reco::TrackCollection> tracks;
00528    iEvent.getByLabel(trackSrc_, tracks);
00529 
00530    // cout << " TRACKS Size " << tracks->size() << endl;
00531 
00532    if(!tracks.isValid()){
00533      cout << "Error! Can't get selectTracks!" << endl;
00534      return ;
00535    }
00536    double trackSin_eta_mid = 0;
00537    double trackCos_eta_mid = 0;
00538    double trackSin_eta_pos = 0;
00539    double trackCos_eta_pos = 0;
00540    double trackSin_eta_neg = 0;
00541    double trackCos_eta_neg = 0;
00542    double trackSin_eta_posgap = 0;
00543    double trackCos_eta_posgap = 0;
00544    double trackSin_eta_neggap = 0;
00545    double trackCos_eta_neggap = 0;
00546    double trackSin_eta = 0;
00547    double trackCos_eta = 0;
00548    double trackSin_even = 0;
00549    double trackCos_even = 0;
00550    double trackSin_odd = 0;
00551    double trackCos_odd = 0;
00552    double trackSin_PosCh = 0;
00553    double trackCos_PosCh = 0;
00554    double trackSin_NegCh = 0;
00555    double trackCos_NegCh = 0;
00556    double face=0; 
00557    double s=0;
00558    
00559    for(reco::TrackCollection::const_iterator j = tracks->begin(); j != tracks->end(); j++){
00560     
00561      track_eta = j->eta();
00562      track_phi = j->phi();
00563      track_pt = j->pt();
00564      track_charge = j->charge();
00565      if ((track_phi<1.0) || (track_phi>1.16)){
00566 
00567      if (fabs(track_eta) < 2.0){
00568        trackSin_eta+=sin(2*track_phi);
00569        trackCos_eta+=cos(2*track_phi);
00570      }
00571      s = tracks->size();
00572       
00573      if((track_charge > 0) && (fabs(track_eta) < 0.75) ){
00574        trackSin_PosCh+=sin(2*track_phi);
00575        trackCos_PosCh+=cos(2*track_phi);
00576      }
00577      if((track_charge < 0) && (fabs(track_eta) < 0.75)){
00578        trackSin_NegCh+=sin(2*track_phi);
00579        trackCos_NegCh+=cos(2*track_phi);
00580      }
00581      int rnd = rand();
00582      face = (double)rnd/(double)RAND_MAX;
00583      
00584      if (fabs(track_eta) < 0.75){
00585      if (face < 0.5){
00586        trackSin_even+=sin(2*track_phi);
00587        trackCos_even+=cos(2*track_phi);
00588      }
00589      
00590      else {
00591        trackSin_odd+=sin(2*track_phi);
00592        trackCos_odd+=cos(2*track_phi);
00593      }
00594      }
00595      if(fabs(track_eta)<0.75){
00596        trackSin_eta_mid+=sin(2*track_phi);
00597        trackCos_eta_mid+=cos(2*track_phi);
00598      }
00599    
00600      if((track_eta >= 0.75) && (track_eta < 2.0)){
00601        trackSin_eta_pos+=sin(2*track_phi);
00602        trackCos_eta_pos+=cos(2*track_phi);
00603      }
00604      if((track_eta <= -0.75) && (track_eta > -2.0)){
00605        trackSin_eta_neg+=sin(2*track_phi);
00606        trackCos_eta_neg+=cos(2*track_phi);
00607      }
00608    
00609      if((track_eta >= 1) && (track_eta < 2.0)){
00610        trackSin_eta_posgap+=sin(2*track_phi);
00611        trackCos_eta_posgap+=cos(2*track_phi);
00612      }
00613      if((track_eta <= -1) && (track_eta > -2.0)){
00614        trackSin_eta_neggap+=sin(2*track_phi);
00615        trackCos_eta_neggap+=cos(2*track_phi);
00616      }
00617 
00618    }
00619    }
00620         trackPsi_eta_mid = 0.5*atan2(trackSin_eta_mid,trackCos_eta_mid);
00621         trackPsi_eta_pos = 0.5*atan2(trackSin_eta_pos,trackCos_eta_pos);
00622         trackPsi_eta_neg = 0.5*atan2(trackSin_eta_neg,trackCos_eta_neg);
00623         trackPsi_eta = 0.5*atan2(trackSin_eta,trackCos_eta);
00624         trackPsi_odd = 0.5*atan2(trackSin_odd,trackCos_odd);
00625         trackPsi_even = 0.5*atan2(trackSin_even,trackCos_even);
00626         trackPsi_PosCh = 0.5*atan2(trackSin_PosCh,trackCos_PosCh);
00627         trackPsi_NegCh = 0.5*atan2(trackSin_NegCh,trackCos_NegCh);
00628         trackPsi_eta_posgap = 0.5*atan2(trackSin_eta_posgap,trackCos_eta_posgap);
00629         trackPsi_eta_neggap = 0.5*atan2(trackSin_eta_neggap,trackCos_eta_neggap);
00630 
00631 //  cout << "  "<< trackPsi_eta_mid << " "<< trackPsi_odd << endl;
00632 
00633    std::auto_ptr<EvtPlaneCollection> evtplaneOutput(new EvtPlaneCollection);
00634     
00635    EvtPlane ecalPlane(ugol2[6],s1[6],s2[6],"Ecal");
00636    EvtPlane hcalPlane(ugol2[3],s1[3],s2[3],"Hcal");
00637    EvtPlane hfPlane(ugol2[11],s1[11],s2[11],"HF");
00638    EvtPlane caloPlane(ugol2[0],s1[0],s2[0],"Calo");
00639    EvtPlane calohfPlane(ugol2[12],s1[12],s2[12],"CaloHF");
00640 
00641    EvtPlane ecalpPlane(ugol2[13],s1[13],s2[13],"EcalP");
00642    EvtPlane ecalmPlane(ugol2[14],s1[14],s2[14],"EcalM");
00643    EvtPlane hcalpPlane(ugol2[15],s1[15],s2[15],"HcalP");
00644    EvtPlane hcalmPlane(ugol2[16],s1[16],s2[16],"HcalM");
00645    EvtPlane hfpPlane(ugol2[9],s1[9],s2[9],"HFp");
00646    EvtPlane hfmPlane(ugol2[10],s1[10],s2[10],"HFm");
00647    EvtPlane calopPlane(ugol2[17],s1[17],s2[17],"CaloP");
00648    EvtPlane calomPlane(ugol2[18],s1[18],s2[18],"CaloM");
00649    EvtPlane calohfpPlane(ugol2[19],s1[19],s2[19],"CaloHFP");
00650    EvtPlane calohfmPlane(ugol2[20],s1[20],s2[20],"CaloHFM");
00651    
00652    EvtPlane hf1Plane(ugol2[23],s1[23],s2[23],"HF1");
00653    EvtPlane hfp1Plane(ugol2[21],s1[21],s2[21],"HFp1");
00654    EvtPlane hfm1Plane(ugol2[22],s1[22],s2[22],"HFm1");
00655 
00656    EvtPlane etecalPlane(etugol2[6],ets1[6],ets2[6],"etEcal");
00657    EvtPlane ethcalPlane(etugol2[3],ets1[3],ets2[3],"etHcal");
00658    EvtPlane ethfPlane(etugol2[11],ets1[11],ets2[11],"etHF");
00659    EvtPlane etcaloPlane(etugol2[0],ets1[0],ets2[0],"etCalo");
00660    EvtPlane etcalohfPlane(etugol2[12],ets1[12],ets2[12],"etCaloHF");
00661 
00662    EvtPlane etecalpPlane(etugol2[13],ets1[13],ets2[13],"etEcalP");
00663    EvtPlane etecalmPlane(etugol2[14],ets1[14],ets2[14],"etEcalM");
00664    EvtPlane ethcalpPlane(etugol2[15],ets1[15],ets2[15],"etHcalP");
00665    EvtPlane ethcalmPlane(etugol2[16],ets1[16],ets2[16],"etHcalM");
00666    EvtPlane ethfpPlane(etugol2[9],ets1[9],ets2[9],"etHFp");
00667    EvtPlane ethfmPlane(etugol2[10],ets1[10],ets2[10],"etHFm");
00668    EvtPlane etcalopPlane(etugol2[17],ets1[17],ets2[17],"etCaloP");
00669    EvtPlane etcalomPlane(etugol2[18],ets1[18],ets2[18],"etCaloM");
00670    EvtPlane etcalohfpPlane(etugol2[19],ets1[19],ets2[19],"etCaloHFP");
00671    EvtPlane etcalohfmPlane(etugol2[20],ets1[20],ets2[20],"etCaloHFM");
00672 
00673    EvtPlane ethf1Plane(etugol2[23],ets1[23],ets2[23],"etHF1");
00674    EvtPlane ethfp1Plane(etugol2[21],ets1[21],ets2[21],"etHFp1");
00675    EvtPlane ethfm1Plane(etugol2[22],ets1[22],ets2[22],"etHFm1");
00676 
00677 
00678    EvtPlane EvtPlaneFromTracksMidEta(trackPsi_eta_mid,trackSin_eta_mid,trackCos_eta_mid,"EvtPlaneFromTracksMidEta");
00679    EvtPlane EvtPlaneFromTracksPosEta(trackPsi_eta_pos,trackSin_eta_pos,trackCos_eta_pos,"EvtPlaneFromTracksPosEta");
00680    EvtPlane EvtPlaneFromTracksNegEta(trackPsi_eta_neg,trackSin_eta_neg,trackCos_eta_neg,"EvtPlaneFromTracksNegEta");
00681    EvtPlane EvtPlaneFromTracksEta(trackPsi_eta,trackSin_eta,trackCos_eta,"EvtPlaneFromTracksEta");  
00682    EvtPlane EvtPlaneFromTracksOdd(trackPsi_odd,trackSin_odd,trackCos_odd,"EvtPlaneFromTracksOdd");
00683    EvtPlane EvtPlaneFromTracksEven(trackPsi_even,trackSin_even,trackCos_even,"EvtPlaneFromTracksEven");
00684    EvtPlane EvtPlaneFromTracksPosCh(trackPsi_PosCh,trackSin_PosCh,trackCos_PosCh,"EvtPlaneFromTracksPosCh");
00685    EvtPlane EvtPlaneFromTracksNegCh(trackPsi_NegCh,trackSin_NegCh,trackCos_NegCh,"EvtPlaneFromTracksNegCh");
00686    EvtPlane EvtPTracksPosEtaGap(trackPsi_eta_posgap,trackSin_eta_posgap,trackCos_eta_posgap,"EvtPTracksPosEtaGap");
00687    EvtPlane EvtPTracksNegEtaGap(trackPsi_eta_neggap,trackSin_eta_neggap,trackCos_eta_neggap,"EvtPTracksNegEtaGap");
00688 
00689    if(useTrackMidEta_) evtplaneOutput->push_back(EvtPlaneFromTracksMidEta);
00690    if(useTrackPosEta_) evtplaneOutput->push_back(EvtPlaneFromTracksPosEta);
00691    if(useTrackNegEta_) evtplaneOutput->push_back(EvtPlaneFromTracksNegEta);
00692    if(useTrackEta_) evtplaneOutput->push_back(EvtPlaneFromTracksEta);
00693    if(useTrackOdd_) evtplaneOutput->push_back(EvtPlaneFromTracksOdd);
00694    if(useTrackEven_) evtplaneOutput->push_back(EvtPlaneFromTracksEven);
00695    if(useTrackPosCh_) evtplaneOutput->push_back(EvtPlaneFromTracksPosCh);
00696    if(useTrackNegCh_) evtplaneOutput->push_back(EvtPlaneFromTracksNegCh);
00697    if(useTrackPosEtaGap_) evtplaneOutput->push_back(EvtPTracksPosEtaGap);
00698    if(useTrackNegEtaGap_) evtplaneOutput->push_back(EvtPTracksNegEtaGap);
00699 
00700  
00701    if(useECAL_) 
00702     {
00703      evtplaneOutput->push_back(ecalPlane);
00704      evtplaneOutput->push_back(ecalpPlane);
00705      evtplaneOutput->push_back(ecalmPlane);
00706      evtplaneOutput->push_back(etecalPlane);
00707      evtplaneOutput->push_back(etecalpPlane);
00708      evtplaneOutput->push_back(etecalmPlane);
00709     }
00710     
00711    if(useHCAL_) 
00712     {
00713      evtplaneOutput->push_back(hcalPlane);
00714      evtplaneOutput->push_back(hcalpPlane);
00715      evtplaneOutput->push_back(hcalmPlane);
00716      evtplaneOutput->push_back(ethcalPlane);
00717      evtplaneOutput->push_back(ethcalpPlane);
00718      evtplaneOutput->push_back(ethcalmPlane);
00719     }
00720     
00721    if(useECAL_ && useHCAL_) 
00722     { 
00723      evtplaneOutput->push_back(caloPlane);
00724      evtplaneOutput->push_back(hfPlane);
00725      evtplaneOutput->push_back(hfpPlane);
00726      evtplaneOutput->push_back(hfmPlane);
00727      evtplaneOutput->push_back(calohfPlane);
00728      evtplaneOutput->push_back(calopPlane);
00729      evtplaneOutput->push_back(calomPlane);
00730      evtplaneOutput->push_back(calohfpPlane);
00731      evtplaneOutput->push_back(calohfmPlane);
00732      evtplaneOutput->push_back(etcaloPlane);
00733      evtplaneOutput->push_back(ethfPlane);
00734      evtplaneOutput->push_back(ethfpPlane);
00735      evtplaneOutput->push_back(ethfmPlane);
00736      evtplaneOutput->push_back(etcalohfPlane);
00737      evtplaneOutput->push_back(etcalopPlane);
00738      evtplaneOutput->push_back(etcalomPlane);
00739      evtplaneOutput->push_back(etcalohfpPlane);
00740      evtplaneOutput->push_back(etcalohfmPlane);    
00741      evtplaneOutput->push_back(hf1Plane);
00742      evtplaneOutput->push_back(hfp1Plane);
00743      evtplaneOutput->push_back(hfm1Plane);
00744      evtplaneOutput->push_back(ethf1Plane);
00745      evtplaneOutput->push_back(ethfp1Plane);
00746      evtplaneOutput->push_back(ethfm1Plane);
00747     
00748      }
00749 
00750    iEvent.put(evtplaneOutput, "recoLevel");
00751 
00752 
00753 // cout << "  "<< planeA << endl;
00754  
00755 }
00756 
00757 // ------------ method called once each job just before starting event loop  ------------
00758 void 
00759 EvtPlaneProducer::beginJob()
00760 {
00761 }
00762 
00763 // ------------ method called once each job just after ending the event loop  ------------
00764 void 
00765 EvtPlaneProducer::endJob() {
00766 }
00767 
00768 //define this as a plug-in
00769 DEFINE_FWK_MODULE(EvtPlaneProducer);
/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/RecoHI/HiEvtPlaneAlgos/src/EvtPlaneProducer.cc
