Changeset 377 in svn for trunk/src/VeryForward.cc

Timestamp:

May 12, 2009, 9:41:00 AM (15 years ago)

Author:

Xavier Rouby

Message:

new PDG table ; switch back to old ZDC and old RomanPots

File:

• trunk/src/VeryForward.cc (modified) (7 diffs)

trunk/src/VeryForward.cc

@@ -32 +32 @@
 #include "VeryForward.h"
 #include "H_RomanPot.h"
+#include "PdgParticle.h"
 #include <iostream>
 #include<cmath>
@@ -94 +95 @@
   beamline1->fill(DET->RP_beam1Card,1,DET->RP_IP_name);                                
   beamline1->offsetElements(DET->RP_offsetEl_s,-DET->RP_offsetEl_x);
-  H_RomanPot * rp220_1 = new H_RomanPot("rp220_1",DET->RP_220_s,DET->RP_220_x*(1E6)); // RP 220m, 2mm, beam 1  
-  H_RomanPot * rp420_1 = new H_RomanPot("rp420_1",DET->RP_420_s,DET->RP_420_x*(1E6)); // RP 420m, 4mm, beam 1  
+  H_RomanPot * rp220_1 = new H_RomanPot("rp220_1",DET->RP_220_s,DET->RP_220_x*1E6); // RP 220m, 2mm, beam 1  
+  H_RomanPot * rp420_1 = new H_RomanPot("rp420_1",DET->RP_420_s,DET->RP_420_x*1E6); // RP 420m, 4mm, beam 1  
   beamline1->add(rp220_1); 
   beamline1->add(rp420_1);
@@ -101 +102 @@
   beamline2->fill(DET->RP_beam2Card,-1,DET->RP_IP_name);                              
   beamline2->offsetElements(DET->RP_offsetEl_s,+DET->RP_offsetEl_x);
-  H_RomanPot * rp220_2 = new H_RomanPot("rp220_2",DET->RP_220_s,DET->RP_220_x*(1E6));// RP 220m, 2mm, beam 2     
-  H_RomanPot * rp420_2 = new H_RomanPot("rp420_2",DET->RP_420_s,DET->RP_420_x*(1E6));// RP 420m, 4mm, beam 2     
+  H_RomanPot * rp220_2 = new H_RomanPot("rp220_2",DET->RP_220_s,DET->RP_220_x*1E6);// RP 220m, 2mm, beam 2     
+  H_RomanPot * rp420_2 = new H_RomanPot("rp420_2",DET->RP_420_s,DET->RP_420_x*1E6);// RP 420m, 4mm, beam 2     
   beamline2->add(rp220_2);
   beamline2->add(rp420_2); 
@@ -110 +111 @@
 
 
+
 void VeryForward::ZDC(ExRootTreeWriter *treeWriter, ExRootTreeBranch *branchZDC, TRootGenParticle *particle)
 {
   TRootZdcHits *elementZdc;
   float energy = particle->E;
-  //TLorentzVector genMomentum; 
 
   // Zero degree calorimeter, for forward neutrons and photons
@@ -121 +122 @@
                            && fabs(particle->Eta) > DET->VFD_min_zdc ) {
     elementZdc = (TRootZdcHits*) branchZDC->NewEntry();
-  
+
+     
+    // for compatibility with 'old' version
+    TLorentzVector genMomentum;
+    genMomentum.SetPxPyPzE(particle->Px, particle->Py, particle->Pz, particle->E);
+    elementZdc->Set(genMomentum); 
+    // ******************
+    
+ 
+    //particle->print();
+ 
     // 1) energy smearing
     float energyS = -1.; 
@@ -158 +169 @@
 
     // 4) object nature : e.m. (photon) or had (neutron) ?
-    elementZdc->hadronic_hit = (bool) (particle->PID==pN);
+    //elementZdc->hadronic_hit = (bool) (particle->PID==pN);
   }
-  
-}
+}
+
+
 void VeryForward::RomanPots(ExRootTreeWriter *treeWriter, ExRootTreeBranch *branchRP220,ExRootTreeBranch *branchFP420,TRootGenParticle *particle) 
 {
-  int pid=particle->PID;
-  
+  float charge = particle->Charge, mass = particle->M;
+  if (mass<-999) { // unitialised!
+          PdgParticle pdg_part = DET->PDGtable[particle->PID];
+          charge  = pdg_part.charge();
+          mass    = pdg_part.mass();
+  }
+  //if(particle->Charge!=1) return; // only particles with Q=+1 can hope to reach RP200/FP420
+ 
   TRootRomanPotHits* elementRP220;
-  TRootForwardTaggerHits* elementFP420;
+  //TRootForwardTaggerHits* elementFP420;
+  TRootRomanPotHits* elementFP420;
   
   TLorentzVector genMomentum;
   genMomentum.SetPxPyPzE(particle->Px, particle->Py, particle->Pz, particle->E);
-  // if forward proton
-  if( (pid == pP)  && (particle->Status == 1) &&  (fabs(genMomentum.Eta()) > DET->CEN_max_calo_fwd) )
+  
+  // to go faster, why not rejecting particles already going into the ZDC?
+  if( (particle->Status == 1) &&  (fabs(genMomentum.Eta()) > DET->CEN_max_calo_fwd) )
     {
-      // !!!!!!!! put here particle->CHARGE and particle->MASS
-      H_BeamParticle p1; /// put here particle->CHARGE and particle->MASS
-      p1.smearAng();
-      p1.smearPos();
+      //cout << "VeryForward :: M = " << mass << "\t Q = " << charge << "\t\t " << particle->PID << endl;
+      H_BeamParticle p1(mass,charge); /// put here particle->CHARGE and particle->MASS
+      p1.smearAng();   p1.smearPos(); // vertex smearing
       p1.setPosition(p1.getX()+DET->RP_cross_x,p1.getY()+DET->RP_cross_y,p1.getTX()-1*kickers_on*DET->RP_cross_ang,p1.getTY(),0);
       p1.set4Momentum(particle->Px,particle->Py,particle->Pz,particle->E);
@@ -222 +241 @@
         } else if (p1.getStoppingElement()->getName()=="rp420_1" || p1.getStoppingElement()->getName()=="rp420_2") {
           p1.propagate(DET->RP_420_s);
-          elementFP420 = (TRootForwardTaggerHits*) branchFP420->NewEntry();
+          //elementFP420 = (TRootForwardTaggerHits*) branchFP420->NewEntry();
+          elementFP420 = (TRootRomanPotHits*) branchFP420->NewEntry();
           elementFP420->X  = (1E-6)*p1.getX();  // [m]
           elementFP420->Y  = (1E-6)*p1.getY();  // [m]