/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 | Private Types | Private Member Functions | Private Attributes
DsPmtSensDet Class Reference

#include <DsPmtSensDet.h>

Collaboration diagram for DsPmtSensDet:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 DsPmtSensDet (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~DsPmtSensDet ()
virtual void Initialize (G4HCofThisEvent *HCE)
virtual void EndOfEvent (G4HCofThisEvent *HCE)
virtual bool ProcessHits (G4Step *step, G4TouchableHistory *history)
virtual StatusCode initialize ()
virtual StatusCode finalize ()

Private Types

typedef std::map< short int,
G4DhHitCollection * > 
LocalHitCache

Private Member Functions

void StoreHit (DayaBay::SimPmtHit *hit, int trackid)
const DetectorElement * SensDetElem (const G4TouchableHistory &hist)
int SensDetId (const DetectorElement &de)
double SensDetQE (G4LogicalVolume *logvol, double energy)

Private Attributes

std::vector< std::string > m_cathodeLogVols
 Properties:
std::vector< std::string > m_sensorStructures
 SensorStructures : names of paths in TDS in which to search for sensor detector elements using this sensitive detector.
std::string m_idParameter
 PackedIdParameterName : name of user paramater of the counted detector element which holds the packed, globally unique PMT ID.
std::string m_t2deName
 TouchableToDetelem : the ITouchableToDetectorElement to use to resolve sensor ID.
ITouchableToDetectorElementm_t2de
double m_qeScale
 QEScale: Upward adjustment of DetSim efficiency to allow PMT-to-PMT efficiency variation in the electronics simulation.
bool m_ConvertWeightToEff
std::string m_qeffParamName
 QEffParameterName : name of user parameter in the photo cathode volume that holds the quantum efficiency tabproperty.
LocalHitCache m_hc
Rndm::Numbers m_uni

Detailed Description

Definition at line 26 of file DsPmtSensDet.h.


Member Typedef Documentation

typedef std::map<short int,G4DhHitCollection*> DsPmtSensDet::LocalHitCache [private]

Definition at line 87 of file DsPmtSensDet.h.


Constructor & Destructor Documentation

DsPmtSensDet::DsPmtSensDet ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

Definition at line 56 of file DsPmtSensDet.cc.

    : G4VSensitiveDetector(name)
    , GiGaSensDetBase(type,name,parent)
    , m_t2de(0)
{
    info() << "DsPmtSensDet (" << type << "/" << name << ") created" << endreq;

    declareProperty("CathodeLogicalVolume",
                    m_cathodeLogVols,
                    "Photo-Cathode logical volume to which this SD is attached.");

    declareProperty("TouchableToDetelem", m_t2deName = "TH2DE",
                    "The ITouchableToDetectorElement to use to resolve sensor.");

    declareProperty("SensorStructures",m_sensorStructures,
                    "TDS Paths in which to look for sensor detector elements"
                    " using this sensitive detector");

    declareProperty("PackedIdPropertyName",m_idParameter="PmtID",
                    "The name of the user property holding the PMT ID.");

    declareProperty("QEffParameterName",m_qeffParamName="EFFICIENCY",
                    "name of user parameter in the photo cathode volume that"
                    " holds the quantum efficiency tabproperty");

    declareProperty("QEScale",m_qeScale=1.0 / 0.9,
                    "Upward scaling of the quantum efficiency by inverse of mean PMT-to-PMT efficiency in electronics simulation.");

    declareProperty("ConvertWeightToEff", m_ConvertWeightToEff=false,            
                    "Treat to the optical photon weight as to preliminary applied QE."
                    "Will affect only the primary photons (GtDiffuserBallTool, etc.).");
    
    m_cathodeLogVols.push_back("/dd/Geometry/PMT/lvPmtHemiCathode");
    m_cathodeLogVols.push_back("/dd/Geometry/PMT/lvHeadonPmtCathode");
}
DsPmtSensDet::~DsPmtSensDet ( ) [virtual]

Definition at line 94 of file DsPmtSensDet.cc.

{
}

Member Function Documentation

void DsPmtSensDet::Initialize ( G4HCofThisEvent *  HCE) [virtual]

Definition at line 195 of file DsPmtSensDet.cc.

{
    m_hc.clear();

    G4DhHitCollection* hc = new G4DhHitCollection(SensitiveDetectorName,collectionName[0]);
    m_hc[0] = hc;
    int hcid = G4SDManager::GetSDMpointer()->GetCollectionID(hc);
    hce->AddHitsCollection(hcid,hc);

    for (int isite=0; site_ids[isite] >= 0; ++isite) {
        for (int idet=0; detector_ids[idet] >= 0; ++idet) {
            DayaBay::Detector det(site_ids[isite],detector_ids[idet]);

            if (det.bogus()) continue;

            string name=det.detName();
            G4DhHitCollection* hc = new G4DhHitCollection(SensitiveDetectorName,name.c_str());
            short int id = det.siteDetPackedData();
            m_hc[id] = hc;

            int hcid = G4SDManager::GetSDMpointer()->GetCollectionID(hc);
            hce->AddHitsCollection(hcid,hc);
            debug() << "Add hit collection with hcid=" << hcid << ", cached ID=" 
                    << (void*)id 
                    << " name= \"" << SensitiveDetectorName << "/" << name << "\"" 
                    << endreq;
        }
    }

    debug() << "DsPmtSensDet Initialize, made "
           << hce->GetNumberOfCollections() << " collections"
           << endreq;
    
}
void DsPmtSensDet::EndOfEvent ( G4HCofThisEvent *  HCE) [virtual]

Definition at line 540 of file DsPmtSensDet.cc.

{
    debug() << "Cache has " << m_hc.size() << " collections" << endreq;

    int ncols = hce->GetNumberOfCollections();
    debug() << "DsPmtSensDet EndOfEvent " << ncols << " collections.";
    //    for (int ind=0; ind<ncols; ++ind) {
    //  G4VHitsCollection* hc = hce->GetHC(ind);
    //  info () << ind << ": " 
    //          << hc->GetSDname() << "//" << hc->GetName() << " has " 
    //          << hc->GetSize() << " hits\n";
    //}
    //info() << endreq;


    int tothits = 0;
    for (int ind=0; ind<ncols; ++ind) {
      G4VHitsCollection* hc = hce->GetHC(ind);
      if ( hc->GetSize() > 0) 
        {
          if ( tothits == 0) debug() << endreq;
          debug() << ind << ": " 
                  << hc->GetSDname() << "//" << hc->GetName() << " has " 
                  << hc->GetSize() << " hits" << endreq;
        }
      tothits += hc->GetSize() ;
    }
    if ( tothits == 0 ) debug() << " No hits found in " << ncols << " collections."  << endreq;
}
bool DsPmtSensDet::ProcessHits ( G4Step *  step,
G4TouchableHistory *  history 
) [virtual]

Definition at line 318 of file DsPmtSensDet.cc.

{
    //if (!step) return false; just crash for now if not defined

    // Find out what detector we are in (ADx, IWS or OWS)
    G4StepPoint* preStepPoint = step->GetPreStepPoint();

    double energyDep = step->GetTotalEnergyDeposit();

    if (energyDep <= 0.0) {
        //debug() << "Hit energy too low: " << energyDep/CLHEP::eV << endreq;
        return false;
    }

    const G4TouchableHistory* hist = 
        dynamic_cast<const G4TouchableHistory*>(preStepPoint->GetTouchable());
    if (!hist or !hist->GetHistoryDepth()) {
        error() << "ProcessHits: step has no or empty touchable history" << endreq;
        return false;
    }

    const DetectorElement* de = this->SensDetElem(*hist);
    if (!de) return false;

    // wangzhe QE calculation starts here.
    int pmtid = this->SensDetId(*de);
    DayaBay::Detector detector(pmtid);
    G4Track* track = step->GetTrack();
    double weight = track->GetWeight();

    // Make sure the energy deposit is caused by an optical photon.
    // In the future we may decide that energy deposits in the photocathode by other particles
    // are capable of generating a PMT signal. In that case, we will need a model for the
    // effects. djaffe

    if ( track->GetDefinition() != G4OpticalPhoton::Definition() ) 
      {
        debug() << "Track definition " << track->GetDefinition()->GetParticleName() 
                << "This track is not an optical photon. abort further processing " << endreq;
        return false;
      }

    G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(track->GetUserInformation());
    DsPhotonTrackInfo* userinfo = composite?dynamic_cast<DsPhotonTrackInfo*>
                                             ( composite->GetPhotonTrackInfo() ):0;
    DsPhotonTrackInfo::QEMode qemode=userinfo?userinfo->GetMode():DsPhotonTrackInfo::kQENone;
    // printf("Got photon %i (w=%.4g)-> ti=%i, mode: %i, pqe=%g, reemitted=%i\n", 
    //          track->GetTrackID(), weight, 
    //          (bool)userinfo, qemode, userinfo?userinfo->GetQE():1., userinfo?userinfo->GetReemitted():0);//gonchar
    debug() << "Got photon " << track->GetTrackID() << " weight "  << weight 
            << " qemode " << qemode << " DsPhotonTrackInfo::kQEWater " << DsPhotonTrackInfo::kQEWater
            << endreq;
    if ( qemode!=DsPhotonTrackInfo::kQEWater ){
        double qe = this->SensDetQE(preStepPoint->GetPhysicalVolume()->GetLogicalVolume(),energyDep);
        debug() << "QE(initial) " << qe << " will be multiple by m_qeScale " << m_qeScale << endreq;
        // Scale total efficiency upward to allow for PMT-to-PMT
        // efficiency variation applied in the electronics simulation.
        // The QEScale used here should be the inverse of the mean
        // PMT efficiency applied in the electronics simulation.
        qe *= m_qeScale;
          
        // For now, hard code "extra" decrease in efficiency for water shield PMTs to match G4dyb.
        if (detector.detectorId() == DetectorId::kIWS || detector.detectorId() == DetectorId::kOWS) {
          qe *= 0.6;
        }
        if (qe < 0.) return false;

        if (userinfo) {
            if (qemode==DsPhotonTrackInfo::kQEPreScale){  
                qe/=userinfo->GetQE(); 
            }        
        }
        else if ( m_ConvertWeightToEff ){
            // to use with DiffuserBallTool
            // since it can not set userinfo->GetQE
            qe*=weight;
            weight=1.;
        }

        if (qe > 1.) {
            warning() << "Non physical QE > 1." << endreq;
        }

        double weight1=weight*qe; 
        debug() << " weight " << weight << " qe " << qe << " weight1 " << weight1 << endreq;
        if (weight1>1.5){
            //reduce the weight effectively according to the qe
            if (weight1 > 10.) {
                weight1 = G4RandGauss::shoot(weight1, sqrt(weight1));
            }
            else {
                weight1 = G4Poisson(weight1);
            }
            if (weight1<=0.) return true;
            weight=weight1;
            debug() << "After applying QE, effective weight is " << weight << endreq;
        } 
        else {
            // apply QE
            double looser = m_uni();
            debug()<<"Apply QE: if qe= "<<qe<<" <=  looser= "<<looser<< " then no PE generated" << endreq;
            if (looser >= qe) {
              //verbose() << "Failed quantum efficiency " << looser << " > " << qe << endreq;
              return true;
            }
            // wz QE calculation and sampling end here. PE is going to be created next.
        }
    }
    else {
        // QE is already applied for Cerenkov in water pool
        // printf("Got photon %i -> ti=%i, mode: %i, pqe=%g, reemitted=%i\n", track->GetTrackID(), (bool)userinfo, qemode, userinfo?userinfo->GetQE():1., userinfo?userinfo->GetReemitted():0);
        // printf("nothing is done\n");
      debug() << " Water pool photon. QE already applied" << endreq;
    }

    double wavelength = CLHEP::twopi*CLHEP::hbarc/energyDep;

#if 1
    if (!pmtid) {
        warning () << "Dumping TouchableHistory:\n";
        for (int ind=0; ind<hist->GetHistoryDepth(); ++ind) {
            warning () << "\t (" << ind << ") " << hist->GetVolume(ind)->GetName() 
                    << " #" << hist->GetReplicaNumber(ind)
                    << " physvol @ " << (void*)hist->GetVolume(ind)
                    << "\n";
        }
        warning() << endreq;
    }

    verbose() << "Dumping TouchableHistory:\n";
    for (int ind=0; ind<hist->GetHistoryDepth(); ++ind) {
        verbose() << "\t (" << ind << ") " << hist->GetVolume(ind)->GetName() 
                  << " #" << hist->GetReplicaNumber(ind)
                  << " physvol @ " << (void*)hist->GetVolume(ind)
                  << "\n";
    }
    verbose() << endreq;
#endif


    DayaBay::SimPmtHit* sphit = new DayaBay::SimPmtHit();

    // base hit

    // Time since event created
    sphit->setHitTime(preStepPoint->GetGlobalTime());

    //#include "G4NavigationHistory.hh"

    const G4AffineTransform& trans = hist->GetHistory()->GetTopTransform();
    const G4ThreeVector& global_pos = preStepPoint->GetPosition();
    G4ThreeVector pos = trans.TransformPoint(global_pos);
    sphit->setLocalPos(pos);
    sphit->setSensDetId(pmtid);
    
    // pmt hit
    // sphit->setDir(...);       // for now
    G4ThreeVector pol = trans.TransformAxis(track->GetPolarization());
    pol = pol.unit();
    G4ThreeVector dir = trans.TransformAxis(track->GetMomentum());
    dir = dir.unit();
    sphit->setPol(pol);
    sphit->setDir(dir);
    sphit->setWavelength(wavelength);
    sphit->setType(0);
    // G4cerr<<"PMT: set hit weight "<<weight<<G4endl; //gonchar
    sphit->setWeight(weight);

    // Dump info on the hit and touchable history
#if 0
    debug() << "hit: id=" << (void*)pmtid << ", "
            << wavelength/CLHEP::nm << " nm, "
            << "pol=[" << pol[0] << "," << pol[1] << "," << pol[2] << "] "
            << "pos=[" << pos[0] << "," << pos[1] << "," << pos[2] << "] "
            << endreq;

    for (int ind=0; ind<hist->GetHistoryDepth(); ++ind) {
        verbose () << "\t (" << ind << ") " << hist->GetVolume(ind)->GetName() 
                << " #" << hist->GetReplicaNumber(ind)
                << " physvol @ " << (void*)hist->GetVolume(ind)
                << "\n";
    }
    verbose() << endreq;
#endif


    int trackid = track->GetTrackID();
    this->StoreHit(sphit,trackid);
    debug() << "Stored photon " << trackid << " weight " << weight << " pmtid " << (void*)pmtid << " wavelength(nm) " << wavelength/CLHEP::nm << endreq;
    return true;
}
StatusCode DsPmtSensDet::initialize ( ) [virtual]

Definition at line 98 of file DsPmtSensDet.cc.

{
    this->GiGaSensDetBase::initialize();
    info() << "DsPmtSensDet initialize" << endreq;

    m_t2de = tool<ITouchableToDetectorElement>(m_t2deName);

    // define collections.

    // first a catch-all
    collectionName.insert("unknown");
    for (int isite=0; site_ids[isite] >= 0; ++isite) {
        Site::Site_t site = site_ids[isite];

        for (int idet=0; detector_ids[idet] >= 0; ++idet) {
            DetectorId::DetectorId_t detid = detector_ids[idet];

            DayaBay::Detector det(site,detid);

            if (det.bogus()) continue;

            string name=det.detName();
            collectionName.insert(name.c_str());
            debug() << "Inserted collectionName " << collectionName.size()-1 
                    << " = \"" << name << "\"" 
                    << "(" << isite << "," << idet << "): "
                    << "(" << site << "," << detid << ")"
                    << endreq;
        }
    }

    // Set up a fast sim model for PMT (see ExN05).  Note, this is
    // kind of circular, but I find no other way to attach a fast sim
    // to a sens det.

    //const G4LogicalVolume* lv = GiGaVolumeUtils::findLVolume(m_cathodeLogVol);
    IGiGaGeomCnvSvc *geoSvc = 0;
    if (service("GiGaGeo",geoSvc,true).isFailure()) {
        fatal() << "Failed to get GiGa Geometry service" << endreq;
        return StatusCode::FAILURE;        
    }
    std::vector<std::string>::iterator volumeIter, 
      volumeEnd = m_cathodeLogVols.end();
    int PMTii = 1;
    char str[10] = "";

    for(volumeIter = m_cathodeLogVols.begin(); volumeIter != volumeEnd; 
        volumeIter++){
      std::string cathodeLogVol = *volumeIter;
      G4LogicalVolume* lv = geoSvc->g4LVolume( cathodeLogVol );
      if (!lv) {
        error() << "Failed to get " << cathodeLogVol << endreq;
        return StatusCode::FAILURE;
      }
      std::string PMT_region_name = "PMT_region";
      std::string PmtModel_name = "PmtModel";
      sprintf(str,"%d",PMTii);
      PMT_region_name += str;
      PmtModel_name += str;
      PMTii++;

      G4Region* pmt_reg = new G4Region(PMT_region_name.c_str());
      pmt_reg->AddRootLogicalVolume(const_cast<G4LogicalVolume*>(lv));
      
      vector<double> cuts; 
      cuts.push_back(1.0*mm);
      cuts.push_back(1.0*mm);
      cuts.push_back(1.0*mm);
      pmt_reg->SetProductionCuts(new G4ProductionCuts());
      pmt_reg->GetProductionCuts()->SetProductionCuts(cuts);
      
      DsPmtModel* pmtModel = new DsPmtModel(PmtModel_name.c_str(),pmt_reg);
      pmtModel = 0;
    }

    // set up random numbers
    IRndmGenSvc *rgs = 0;
    if (service("RndmGenSvc",rgs,true).isFailure()) {
        fatal() << "Failed to get random service" << endreq;
        return StatusCode::FAILURE;        
    }
    StatusCode sc = m_uni.initialize(rgs, Rndm::Flat(0,1));
    if (sc.isFailure()) {
        fatal() << "Failed to initialize uniform random numbers" << endreq;
        return StatusCode::FAILURE;
    }

    return StatusCode::SUCCESS;
}
StatusCode DsPmtSensDet::finalize ( ) [virtual]

Definition at line 188 of file DsPmtSensDet.cc.

{
    info() << "DsPmtSensDet finalize" << endreq;
    return this->GiGaSensDetBase::finalize();
}
void DsPmtSensDet::StoreHit ( DayaBay::SimPmtHit hit,
int  trackid 
) [private]

Definition at line 511 of file DsPmtSensDet.cc.

{
    int did = hit->sensDetId();
    DayaBay::Detector det(did);
    short int sdid = det.siteDetPackedData();

    G4DhHitCollection* hc = m_hc[sdid];
    if (!hc) {
        warning() << "Got hit with no hit collection.  ID = " << (void*)did
                  << " which is detector: \"" << DayaBay::Detector(did).detName()
                  << "\". Storing to the " << collectionName[0] << " collection"
                  << endreq;
        sdid = 0;
        hc = m_hc[sdid];
    }

#if 1
    verbose() << "Storing hit PMT: " << (void*)did 
              << " from " << DayaBay::Detector(did).detName()
              << " in hc #"<<  sdid << " = "
              << hit->hitTime()/CLHEP::ns << "[ns] "
              << hit->localPos()/CLHEP::cm << "[cm] "
              << hit->wavelength()/CLHEP::nm << "[nm]"
              << endreq;
#endif

    hc->insert(new G4DhHit(hit,trackid));
}
const DetectorElement * DsPmtSensDet::SensDetElem ( const G4TouchableHistory &  hist) [private]

Definition at line 231 of file DsPmtSensDet.cc.

{
    const IDetectorElement* idetelem = 0;
    int steps=0;

    if (!hist.GetHistoryDepth()) {
        error() << "DsPmtSensDet::SensDetElem given empty touchable history" << endreq;
        return 0;
    }

    StatusCode sc = 
        m_t2de->GetBestDetectorElement(&hist,m_sensorStructures,idetelem,steps);
    if (sc.isFailure()) {      // verbose warning
        warning() << "Failed to find detector element in:\n";
        for (size_t ind=0; ind<m_sensorStructures.size(); ++ind) {
            warning() << "\t\t" << m_sensorStructures[ind] << "\n";
        }
        warning() << "\tfor touchable history:\n";
        for (int ind=0; ind < hist.GetHistoryDepth(); ++ind) {
            warning() << "\t (" << ind << ") " 
                      << hist.GetVolume(ind)->GetName() << "\n";
        }
        warning() << endreq;
        return 0;
    }
    
    return dynamic_cast<const DetectorElement*>(idetelem);
}
int DsPmtSensDet::SensDetId ( const DetectorElement &  de) [private]

Definition at line 260 of file DsPmtSensDet.cc.

{
    const DetectorElement* detelem = &de;

    while (detelem) {
        if (detelem->params()->exists(m_idParameter)) {
            break;
        }
        detelem = dynamic_cast<const DetectorElement*>(detelem->parentIDetectorElement());
    }
    if (!detelem) {
        warning() << "Could not get PMT detector element starting from " << de << endreq;
        return 0;
    }

    return detelem->params()->param<int>(m_idParameter);
}
double DsPmtSensDet::SensDetQE ( G4LogicalVolume *  logvol,
double  energy 
) [private]

Definition at line 277 of file DsPmtSensDet.cc.

{
    G4Material* mat = logvol->GetMaterial();
    if (!mat) {
        warning () << "No material for " << logvol->GetName() << endreq;
        return -1;
    }
    

     G4MaterialPropertiesTable* mattab = mat->GetMaterialPropertiesTable();
    if (mattab) {
        G4MaterialPropertyVector* qevec = mattab->GetProperty(m_qeffParamName.c_str());
        if (qevec) {
        
          verbose() << m_qeffParamName << ":("
                  << qevec->GetMinPhotonEnergy()/CLHEP::eV << " eV," 
                  << qevec->GetMinProperty() << ")-->("
                  << qevec->GetMaxPhotonEnergy()/CLHEP::eV << " eV," 
                  << qevec->GetMaxProperty() << ")"
                  << " particle energy is " << energy/CLHEP::eV
                  << endreq;

          return qevec->GetProperty(energy);

        }
    }
    else {
        debug () << "No material properties in " << logvol->GetName() << endreq;
    }

    int ndaught = logvol->GetNoDaughters();
    for (int ind=0; ind < ndaught; ++ind) {
        G4VPhysicalVolume* physvol = logvol->GetDaughter(ind);
        double qe = this->SensDetQE(physvol->GetLogicalVolume(),energy);
        if (qe < 0) return qe;
    }
    warning() << "All attempts failed to find " << m_qeffParamName 
              << " in " << logvol->GetName() << endreq;
    return -1;
}

Member Data Documentation

std::vector<std::string> DsPmtSensDet::m_cathodeLogVols [private]

Properties:

CathodeLogicalVolumes : name of logical volumes in which this sensitive detector is operating.

Definition at line 48 of file DsPmtSensDet.h.

std::vector<std::string> DsPmtSensDet::m_sensorStructures [private]

SensorStructures : names of paths in TDS in which to search for sensor detector elements using this sensitive detector.

Definition at line 52 of file DsPmtSensDet.h.

std::string DsPmtSensDet::m_idParameter [private]

PackedIdParameterName : name of user paramater of the counted detector element which holds the packed, globally unique PMT ID.

Definition at line 57 of file DsPmtSensDet.h.

std::string DsPmtSensDet::m_t2deName [private]

TouchableToDetelem : the ITouchableToDetectorElement to use to resolve sensor ID.

Definition at line 61 of file DsPmtSensDet.h.

Definition at line 62 of file DsPmtSensDet.h.

double DsPmtSensDet::m_qeScale [private]

QEScale: Upward adjustment of DetSim efficiency to allow PMT-to-PMT efficiency variation in the electronics simulation.

The value should be the inverse of the mean PMT efficiency applied in ElecSim.

Definition at line 68 of file DsPmtSensDet.h.

Definition at line 71 of file DsPmtSensDet.h.

std::string DsPmtSensDet::m_qeffParamName [private]

QEffParameterName : name of user parameter in the photo cathode volume that holds the quantum efficiency tabproperty.

Definition at line 75 of file DsPmtSensDet.h.

Definition at line 88 of file DsPmtSensDet.h.

Rndm::Numbers DsPmtSensDet::m_uni [private]

Definition at line 90 of file DsPmtSensDet.h.


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

Generated on Fri May 16 2014 10:17:58 for DetSim by doxygen 1.7.4