Commit 213b3848 authored by Simon Spannagel's avatar Simon Spannagel
Browse files

Deposition modules: also store total charge generated by MCParticle

parent a212ffa1
Loading
Loading
Loading
Loading
+4 −0
Original line number Diff line number Diff line
@@ -108,6 +108,7 @@ G4bool SensitiveDetectorActionG4::ProcessHits(G4Step* step, G4TouchableHistory*)
    // Update current end point with the current last step
    auto end_position = detector_->getLocalPosition(static_cast<ROOT::Math::XYZPoint>(postStep->GetPosition()));
    track_end_[trackID] = end_position;
    track_charge_[trackID] += charge;

    // Add new deposit if the charge is more than zero
    if(charge == 0) {
@@ -152,6 +153,7 @@ void SensitiveDetectorActionG4::clearEventInfo() {
    track_end_.clear();
    track_pdg_.clear();
    track_time_.clear();
    track_charge_.clear();

    deposit_position_.clear();
    deposit_charge_.clear();
@@ -178,6 +180,7 @@ void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* mess
        ROOT::Math::XYZPoint end_point;
        auto local_end = track_end_.at(track_id);
        auto pdg_code = track_pdg_.at(track_id);
        auto charge = track_charge_.at(track_id);
        auto track_time_global = track_time_.at(track_id);
        auto track_time_local = track_time_global - time_reference;

@@ -185,6 +188,7 @@ void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* mess
        auto global_end = detector_->getGlobalPosition(local_end);
        mc_particles.emplace_back(
            local_begin, global_begin, local_end, global_end, pdg_code, track_time_local, track_time_global);
        mc_particles.back().setTotalDepositedCharge(charge);
        mc_particles.back().setTrack(track_info_manager_->findMCTrack(track_id));
        id_to_particle_[track_id] = mc_particles.size() - 1;

+2 −0
Original line number Diff line number Diff line
@@ -140,6 +140,8 @@ namespace allpix {
        std::map<int, int> track_pdg_;
        // Arrival timestamp of the tracks
        std::map<int, double> track_time_;
        // Total charge by track
        std::map<int, unsigned int> track_charge_;

        // Map from deposit index to track id
        std::vector<int> deposit_to_id_;
+5 −1
Original line number Diff line number Diff line
@@ -176,6 +176,7 @@ void DepositionPointChargeModule::DepositPoint(Event* event, const ROOT::Math::X
    mcparticles.emplace_back(position, position_global, position, position_global, -1, 0., 0.);
    LOG(DEBUG) << "Generated MCParticle at global position " << Units::display(position_global, {"um", "mm"})
               << " in detector " << detector_->getName();
    mcparticles.back().setTotalDepositedCharge(carriers_);

    charges.emplace_back(position, position_global, CarrierType::ELECTRON, carriers_, 0., 0., &(mcparticles.back()));
    charges.emplace_back(position, position_global, CarrierType::HOLE, carriers_, 0., 0., &(mcparticles.back()));
@@ -209,14 +210,17 @@ void DepositionPointChargeModule::DepositLine(Event* event, const ROOT::Math::XY
    auto start_global = detector_->getGlobalPosition(start_local);
    auto end_global = detector_->getGlobalPosition(end_local);

    // Total number of carriers will be:
    auto charge = carriers_ * (end_local.z() - start_local.z()) / step_size_z_;
    // Create MCParticle:
    mcparticles.emplace_back(start_local, start_global, end_local, end_global, -1, 0., 0.);
    LOG(DEBUG) << "Generated MCParticle with start " << Units::display(start_global, {"um", "mm"}) << " and end "
               << Units::display(end_global, {"um", "mm"}) << " in detector " << detector_->getName();
    mcparticles.back().setTotalDepositedCharge(charge);

    // Deposit the charge carriers:
    auto position_local = start_local;
    while(position_local.z() < model->getSensorSize().z() / 2.0) {
    while(position_local.z() < end_local.z()) {
        position_local += ROOT::Math::XYZVector(0, 0, step_size_z_);
        auto position_global = detector_->getGlobalPosition(position_local);

+4 −0
Original line number Diff line number Diff line
@@ -220,6 +220,7 @@ void DepositionReaderModule::run(Event* event) {
    std::map<std::shared_ptr<Detector>, std::vector<int>> mc_particle_code;
    std::map<std::shared_ptr<Detector>, std::vector<double>> mc_particle_time;
    std::map<std::shared_ptr<Detector>, std::vector<int>> mc_particle_parent;
    std::map<std::shared_ptr<Detector>, std::vector<unsigned int>> mc_particle_charge;

    std::map<std::shared_ptr<Detector>, std::vector<int>> particles_to_deposits;
    std::map<std::shared_ptr<Detector>, std::map<int, size_t>> track_id_to_mcparticle;
@@ -307,10 +308,12 @@ void DepositionReaderModule::run(Event* event) {
            mc_particle_time[detector].push_back(time);
            mc_particle_code[detector].push_back(pdg_code);
            mc_particle_parent[detector].push_back(parent_id);
            mc_particle_charge[detector].push_back(charge);
            track_id_to_mcparticle[detector][track_id] = (mc_particle_start[detector].size() - 1);
        } else {
            LOG(DEBUG) << "Found MCParticle with track id " << track_id << ", updating position";
            mc_particle_end[detector].at(iter->second) = global_position;
            mc_particle_charge[detector].at(iter->second) += charge;
        }

        particles_to_deposits[detector].push_back(track_id);
@@ -344,6 +347,7 @@ void DepositionReaderModule::run(Event* event) {

            mc_particles.emplace_back(
                start_local, start_global, end_local, end_global, pdg_code, time - time_reference, time);
            mc_particles.back().setTotalDepositedCharge(mc_particle_charge[detector].at(i));
        }

        for(size_t i = 0; i < mc_particle_size; i++) {