/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Public Attributes
chkIBD15::plotIbdBasics Class Reference
Inheritance diagram for chkIBD15::plotIbdBasics:
Inheritance graph
[legend]
Collaboration diagram for chkIBD15::plotIbdBasics:
Collaboration graph
[legend]

List of all members.

Public Member Functions

def __init__
def initialize
def execute
def finalize

Public Attributes

 hist
 SpillIn_onGd
 SpillOut_onGd
 SpillIn_onH
 SpillOut_onH
 target_de_name
 gds_de_name
 lso_de_name
 coorSvc
 histSvc

Detailed Description

Definition at line 20 of file __init__.py.


Constructor & Destructor Documentation

def chkIBD15::plotIbdBasics::__init__ (   self,
  name 
)

Definition at line 22 of file __init__.py.

00023                            :
00024         GaudiAlgo.__init__(self,name)
00025         self.hist = {}
00026         # Counters
00027         self.SpillIn_onGd = 0
00028         self.SpillOut_onGd = 0
00029         self.SpillIn_onH = 0
00030         self.SpillOut_onH = 0
00031 
00032         return
    

Member Function Documentation

def chkIBD15::plotIbdBasics::initialize (   self)

Definition at line 33 of file __init__.py.

00034                         :
00035         print "Initializing the IBD basic ploter", self.name()
00036         status = GaudiAlgo.initialize(self)
00037         if status.isFailure(): return status
00038         self.target_de_name = '/dd/Structure/AD/db-ade1/db-sst1/db-oil1'
00039         self.gds_de_name = '/dd/Structure/AD/db-gds1'
00040         self.lso_de_name = '/dd/Structure/AD/db-lso1'
00041 
00042         # What services do you need?
00043         self.coorSvc = self.svc('ICoordSysSvc', 'CoordSysSvc')
00044         if not self.coorSvc:
00045             print 'Failed to get CoordSysSvc'
00046             return FAILURE
00047 
00048         self.histSvc = self.svc('ITHistSvc', 'THistSvc')
00049 
00050         self.hist["genRZ"] = TH2D("genRZ", "Generation Vertex R-Z", \
00051                                       100, 0.0, 6.190144, 100, -2.48, 2.48)
00052         status = self.histSvc.regHist('/file1/basics/genRZ', \
00053                                           self.hist["genRZ"])
00054         if status.isFailure(): return status
00055         
00056         self.hist["genXY"] = TH2D("genXY", "Generation Vertex X-Y", \
00057                                       100, -2.48, 2.48, 100, -2.48, 2.48)
00058         status = self.histSvc.regHist('/file1/basics/genXY', \
00059                                           self.hist["genXY"])
00060         if status.isFailure(): return status
00061         
00062         self.hist["HitTime"] = TH1F("HitTime", "Hit Time [ms]",
00063                                     500, 0.0, 2000)
00064         status = self.histSvc.regHist('/file1/basics/HitTime', 
00065                                       self.hist["HitTime"])
00066         if status.isFailure(): return status
00067 
00068         self.hist["pe_E"] = TH2F("pe_E", "pe vs visible E (e+ within GdLS)",
00069                                  100, 0, 10, 700, 0, 1400)
00070         status = self.histSvc.regHist('/file1/basics/pe_E', 
00071                                       self.hist["pe_E"])
00072         if status.isFailure(): return status
00073 
00074         self.hist["pe_E_inLS"] = TH2F("pe_E_inLS", \
00075                                         "pe vs visible E (e+ within LS)", \
00076                                         100, 0, 10, 700, 0, 1400)
00077         status = self.histSvc.regHist('/file1/basics/pe_E_inLS', 
00078                                       self.hist["pe_E_inLS"])
00079         if status.isFailure(): return status
00080 
00081         self.hist["pe_E_LS"] = TH2F("pe_E_LS", \
00082                                         "pe vs visible E (e+ in LS)", \
00083                                         100, 0, 10, 700, 0, 1400)
00084         status = self.histSvc.regHist('/file1/basics/pe_E_LS', 
00085                                       self.hist["pe_E_LS"])
00086         if status.isFailure(): return status
00087 
00088         self.hist["pe_E_all"] = TH2F("pe_E_all", "pe vs visible E (all e+)",
00089                                      100, 0, 10, 700, 0, 1400)
00090         status = self.histSvc.regHist('/file1/basics/pe_E_all', 
00091                                       self.hist["pe_E_all"])
00092         if status.isFailure(): return status
00093         
00094         self.hist["nHits_LS"] = TH1F("nHits_LS", \
00095                                          "nCap Hits in LS", \
00096                                          100, 0, 1400)
00097         status = self.histSvc.regHist('/file1/basics/nHits_LS', 
00098                                       self.hist["nHits_LS"])
00099         if status.isFailure(): return status
00100         
00101         self.hist["nHits_MO"] = TH1F("nHits_MO", \
00102                                          "nCap Hits in MO", \
00103                                          100, 0, 1400)
00104         status = self.histSvc.regHist('/file1/basics/nHits_MO', 
00105                                       self.hist["nHits_MO"])
00106         if status.isFailure(): return status
00107 
00108         self.hist["nHits"] = TH1F("nHits", "Neutron Capture Hits",
00109                                   100, 0, 1400)
00110         status = self.histSvc.regHist('/file1/basics/nHits', 
00111                                       self.hist["nHits"])
00112         if status.isFailure(): return status
00113 
00114         self.hist["nHits_all"] = TH1F("nHits_all", "Neutron Capture Hits",
00115                                   100, 0, 1400)
00116         status = self.histSvc.regHist('/file1/basics/nHits_all', 
00117                                       self.hist["nHits_all"])
00118         if status.isFailure(): return status
00119 
00120         self.hist["nCapPos"] = TH1F("nCapPos", "Neutron Capture Position",
00121                                     310, 0, 6.190144)
00122         status = self.histSvc.regHist('/file1/basics/nCapPos',
00123                                       self.hist["nCapPos"])
00124         if status.isFailure(): return status
00125 
00126         self.hist["eDepInGdLS"] = TH1F("eDepInGdLS", "Deposited Energy [MeV]",
00127                                        20, 0, 20)
00128         status = self.histSvc.regHist('/file1/basics/eDepInGdLS', 
00129                                       self.hist["eDepInGdLS"])
00130         if status.isFailure(): return status
00131 
00132         self.hist["eDepInLS"] = TH1F("eDepInLS", "Deposited Energy [MeV]",
00133                                        20, 0, 20)
00134         status = self.histSvc.regHist('/file1/basics/eDepInLS', 
00135                                       self.hist["eDepInLS"])
00136         if status.isFailure(): return status
00137 
00138         self.hist["drift_GdLS"] = TH1F("drift_GdLS",
00139                                       "Neutron Drift Distance in GdLS [cm]",
00140                                        200, 0, 50)
00141         status = self.histSvc.regHist('/file1/basics/drift_GdLS', 
00142                                       self.hist["drift_GdLS"])
00143         if status.isFailure(): return status
00144 
00145         self.hist["drift_LS"] = TH1F("drift_LS",
00146                                       "Neutron Drift Distance in LS [cm]",
00147                                        200, 0, 50)
00148         status = self.histSvc.regHist('/file1/basics/drift_LS', 
00149                                       self.hist["drift_LS"])
00150         if status.isFailure(): return status
00151         
00152         self.hist["time_GdLS"] = TH1F("time_GdLS",
00153                                       "Neutron Capture time in GdLS [ms]",
00154                                        400, 0, 400)
00155         status = self.histSvc.regHist('/file1/basics/time_GdLS', 
00156                                       self.hist["time_GdLS"])
00157         if status.isFailure(): return status
00158 
00159         self.hist["time_LS"] = TH1F("time_LS",
00160                                       "Neutron Capture Time in LS [ms]",
00161                                        400, 0, 2000)
00162         status = self.histSvc.regHist('/file1/basics/time_LS', 
00163                                       self.hist["time_LS"])
00164         if status.isFailure(): return status
00165 
00166         return SUCCESS

def chkIBD15::plotIbdBasics::execute (   self)

Definition at line 167 of file __init__.py.

00168                      :
00169         print "Executing plotIbdBasics", self.name()
00170         evt = self.evtSvc()
00171         simhdr = evt['/Event/Sim/SimHeader']
00172 #        print "SimHeader: ", simhdr
00173         if simhdr == None:
00174             print "No SimHeader in this ReadOut"
00175             return SUCCESS
00176 
00177         det = self.detSvc(self.target_de_name)
00178         det_gds = self.detSvc(self.gds_de_name)
00179         det_lso = self.detSvc(self.lso_de_name)
00180 
00181         # Unobservables
00182         statshdr = simhdr.unobservableStatistics()
00183         stats = statshdr.stats()
00184 
00185         PID_trk1 = stats["pdgId_Trk1"].sum()
00186         PID_trk2 = stats["pdgId_Trk2"].sum()
00187 
00188         if PID_trk1 != -11 or PID_trk2 != 2112:
00189             print "PID of track 1 is", PID_trk1
00190             print "PID of track 2 is", PID_trk2
00191             print "Not an IBD event."
00192             return SUCCESS
00193 
00194         tGen = stats["t_Trk2"].sum()
00195         xGen = stats["x_Trk2"].sum()
00196         yGen = stats["y_Trk2"].sum()
00197         zGen = stats["z_Trk2"].sum()
00198         
00199         tCap = stats["tEnd_Trk2"].sum()
00200         xCap = stats["xEnd_Trk2"].sum()
00201         yCap = stats["yEnd_Trk2"].sum()
00202         zCap = stats["zEnd_Trk2"].sum()
00203         
00204         # Get underlying DE object
00205         de = self.getDet(self.target_de_name)
00206         de_lso = self.getDet(self.lso_de_name)
00207         de_gds = self.getDet(self.gds_de_name)
00208         if not de:
00209             print 'Failed to get DE',self.target_de_name
00210             return FAILURE
00211         
00212         if not de_lso:
00213             print 'Failed to get DE',self.lso_de_name
00214             return FAILURE
00215         
00216         if not de_gds:
00217             print 'Failed to get DE',self.gds_de_name
00218             return FAILURE
00219         
00220         # Get the AD coordinates of the vertexes
00221         import PyCintex
00222         Gaudi = PyCintex.makeNamespace('Gaudi')
00223         genGlbPoint = Gaudi.XYZPoint(xGen, yGen, zGen)
00224         capGlbPoint = Gaudi.XYZPoint(xCap, yCap, zCap)
00225 #        point = de.geometry().toGlobal(point)
00226         genLclPoint = de.geometry().toLocal(genGlbPoint)
00227         capLclPoint = de.geometry().toLocal(capGlbPoint)
00228         genLclPointLso = de_lso.geometry().toLocal(genGlbPoint)
00229         capLclPointLso = de_lso.geometry().toLocal(capGlbPoint)
00230         genLclPointGds = de_gds.geometry().toLocal(genGlbPoint)
00231         capLclPointGds = de_gds.geometry().toLocal(capGlbPoint)
00232 #        print 'In global coordinate [',gpoint.x(),gpoint.y(),gpoint.z(),']'
00233 
00234         ndrift = ROOT.TVector3(xCap-xGen, yCap-yGen, zCap-zGen)
00235         
00236         capTime = tCap - tGen
00237         capDis = ndrift.Mag()
00238 
00239         R2 = genLclPoint.x()/units.meter * genLclPoint.x()/units.meter + \
00240             genLclPoint.y()/units.meter * genLclPoint.y()/units.meter
00241         
00242         self.hist["genRZ"].Fill(R2, genLclPoint.z()/units.meter)
00243 
00244         self.hist["genXY"].Fill(genLclPoint.x()/units.meter,genLclPoint.y()/units.meter)
00245 
00246         # Find the interesting volumes
00247         genDE = self.coorSvc.coordSysDE(genGlbPoint)
00248         capDE = self.coorSvc.coordSysDE(capGlbPoint)
00249         if not genDE:
00250             print 'Failed to find coordinate system DE for generation'\
00251                 '[',genGlbPoint.x(),genGlbPoint.y(),genGlbPoint.z(),']'
00252             print 'Local: [',genLclPoint.x(),genLclPoint.y(),genLclPoint.z(),']'
00253             return FAILURE
00254         else:
00255             gendmvol = genDE.geometry().belongsToPath(genGlbPoint,-1)
00256 
00257         if not capDE:
00258             print 'Failed to find coordinate system DE for capture'\
00259                 '[',capGlbPoint.x(),capGlbPoint.y(),capGlbPoint.z(),']'
00260             return FAILURE
00261         else:
00262             capdmvol = capDE.geometry().belongsToPath(capGlbPoint,-1)
00263 
00264         import re
00265         genDM = re.split('/', gendmvol).pop()
00266         capDM = re.split('/', capdmvol).pop()
00267         print "Generated in ", genDM
00268         print "Captured in ", capDM
00269 
00270         positronHits = 0
00271         neutronHits = 0
00272         positronTimeCut = 500.
00273         simhits = simhdr.hits()
00274         for detector in simhits.hitDetectors():
00275             hitCollection = simhits.hitsByDetector(detector)
00276             if hitCollection == None:
00277                 print "No hits in ", detector
00278             hits = hitCollection.collection()
00279             for hit in hits:
00280 #                print " PMT", hit.sensDetId(), "hit @ time: ", hit.hitTime()
00281                 self.hist["HitTime"].Fill(hit.hitTime()/units.microsecond)
00282                 if hit.hitTime()/units.nanosecond<positronTimeCut and \
00283                         hit.hitTime()/units.nanosecond < \
00284                         capTime/units.nanosecond:
00285                     positronHits += 1
00286                 else:
00287                     neutronHits += 1
00288 
00289         # Visible energy
00290         vis_trk1 = 1.022 + stats['ke_Trk1'].sum()/units.MeV
00291         self.hist["pe_E_all"].Fill(vis_trk1, positronHits)
00292         if genDM ==  'db-gds1':
00293             self.hist["pe_E"].Fill(vis_trk1, positronHits)
00294         if genDM ==  'db-lso1':
00295             self.hist["pe_E_LS"].Fill(vis_trk1, positronHits)
00296         if re.search('db-lso1',gendmvol):
00297             self.hist["pe_E_inLS"].Fill(vis_trk1, positronHits)
00298 
00299         capTarget = stats["capTarget"].sum()
00300         print "The capture target is ", capTarget
00301 
00302         if capTarget == 1 and \
00303                 capLclPoint.z()/units.m < 1. and \
00304                 capLclPoint.z()/units.m > -1.:
00305             self.hist["nCapPos"].Fill(R2)
00306 
00307         self.hist["nHits_all"].Fill(neutronHits)
00308 
00309         if capDM == 'db-lso1':
00310             self.hist["nHits_LS"].Fill(neutronHits)
00311 
00312         if capDM == 'db-oil1':
00313             self.hist["nHits_MO"].Fill(neutronHits)
00314 
00315         if genDM == 'db-gds1' and capDM == 'db-gds1':
00316             self.hist["nHits"].Fill(neutronHits)
00317             self.hist["time_GdLS"].Fill(capTime/units.microsecond)
00318             self.hist["drift_GdLS"].Fill(capDis/units.cm)
00319 
00320         if genDM == 'db-lso1' and capDM == 'db-lso1':
00321             self.hist["time_LS"].Fill(capTime/units.microsecond)
00322             self.hist["drift_LS"].Fill(capDis/units.cm)
00323 
00324         if genDM == 'db-lso1' and capDM == 'db-oil1':
00325             self.SpillOut_onH += 1
00326 
00327         if genDM == 'db-oil1' and capDM == 'db-lso1':
00328             self.SpillIn_onH += 1
00329 
00330         if genDM == 'db-gds1' and capDM == 'db-lso1':
00331             self.SpillOut_onGd += 1
00332 
00333         if genDM == 'db-lso1' and capDM == 'db-gds1':
00334             self.SpillIn_onGd += 1
00335 
00336         eDepInGdLS = stats["EDepInGdLS"].sum()
00337         self.hist["eDepInGdLS"].Fill(eDepInGdLS/units.MeV)
00338         eDepInLS = stats["EDepInLS"].sum()
00339         self.hist["eDepInLS"].Fill(eDepInLS/units.MeV)
00340         
00341         return SUCCESS
    
def chkIBD15::plotIbdBasics::finalize (   self)

Definition at line 342 of file __init__.py.

00343                       :
00344 
00345         print "Spill-In of on Gd n: ", self.SpillIn_onGd
00346         print "Spill-Out of n from Gd-LS: ", self.SpillOut_onGd
00347         print "Spill-In of n from OIL to LSO: ", self.SpillIn_onH
00348         print "Spill-Out of n from LSO to OIL: ", self.SpillOut_onH
00349 
00350         print "Finalizing ", self.name()
00351         status = GaudiAlgo.finalize(self)
00352         return status


Member Data Documentation

Definition at line 22 of file __init__.py.

Definition at line 22 of file __init__.py.

Definition at line 22 of file __init__.py.

Definition at line 22 of file __init__.py.

Definition at line 22 of file __init__.py.

Definition at line 33 of file __init__.py.

Definition at line 33 of file __init__.py.

Definition at line 33 of file __init__.py.

Definition at line 33 of file __init__.py.

Definition at line 33 of file __init__.py.


The documentation for this class was generated from the following file:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 09:55:05 for MDC09a by doxygen 1.7.4