Commit 14ce5111 authored by Paul Schütze's avatar Paul Schütze
Browse files

DepositionG4: add histogram diplaying the deposited energy per event

parent 5b88f1bb
Loading
Loading
Loading
Loading
+15 −2
Original line number Diff line number Diff line
@@ -330,6 +330,9 @@ void DepositionGeant4Module::run(Event* event) {
        if(output_plots_) {
            double charge = static_cast<double>(Units::convert(sensor->getDepositedCharge(), "ke"));
            charge_per_event_[sensor->getName()]->Fill(charge);

            double deposited_energy = static_cast<double>(Units::convert(sensor->getDepositedEnergy(), "keV"));
            energy_per_event_[sensor->getName()]->Fill(deposited_energy);
        }
    }

@@ -343,6 +346,9 @@ void DepositionGeant4Module::finalize() {
        for(auto& histogram : charge_per_event_) {
            histogram.second->Write();
        }
        for(auto& histogram : energy_per_event_) {
            histogram.second->Write();
        }
    }

    // Print summary or warns if module did not output any charges
@@ -423,14 +429,21 @@ void DepositionGeant4Module::construct_sensitive_detectors_and_fields(double fan
            int nbins = 5 * maximum;

            // Create histograms if needed
            std::string plot_name = "deposited_charge_" + sensitive_detector_action->getName();

            {
                std::lock_guard<std::mutex> lock(histogram_mutex_);
                std::string plot_name = "deposited_charge_" + sensitive_detector_action->getName();

                if(charge_per_event_.find(sensitive_detector_action->getName()) == charge_per_event_.end()) {
                    charge_per_event_[sensitive_detector_action->getName()] = CreateHistogram<TH1D>(
                        plot_name.c_str(), "deposited charge per event;deposited charge [ke];events", nbins, 0, maximum);
                }

                plot_name = "deposited_energy_" + sensitive_detector_action->getName();

                if(energy_per_event_.find(sensitive_detector_action->getName()) == energy_per_event_.end()) {
                    energy_per_event_[sensitive_detector_action->getName()] = CreateHistogram<TH1D>(
                        plot_name.c_str(), "deposited energy per event;deposited energy [keV];events", 200, 0, 1000);
                }
            }
        }
    }
+1 −0
Original line number Diff line number Diff line
@@ -134,6 +134,7 @@ namespace allpix {

        // Vector of histogram pointers for debugging plots
        std::map<std::string, Histogram<TH1D>> charge_per_event_;
        std::map<std::string, Histogram<TH1D>> energy_per_event_;

        // Total deposited charges
        std::atomic_uint total_charges_{0};
+16 −0
Original line number Diff line number Diff line
@@ -129,6 +129,7 @@ G4bool SensitiveDetectorActionG4::ProcessHits(G4Step* step, G4TouchableHistory*)
    // Store relevant quantities to create charge deposits:
    deposit_position_.push_back(deposit_position);
    deposit_charge_.push_back(charge);
    deposit_energy_.push_back(edep);
    deposit_time_.push_back(step_time);
    deposit_to_id_.push_back(trackID);

@@ -152,6 +153,14 @@ unsigned int SensitiveDetectorActionG4::getDepositedCharge() const {
    return deposited_charge_;
}

double SensitiveDetectorActionG4::getTotalDepositedEnergy() const {
    return total_deposited_energy_;
}

double SensitiveDetectorActionG4::getDepositedEnergy() const {
    return deposited_energy_;
}

void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* messenger, Event* event) {

    auto time_reference = std::min_element(track_time_.begin(), track_time_.end(), [](const auto& l, const auto& r) {
@@ -210,6 +219,7 @@ void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* mess

    // Send a deposit message if we have any deposits
    unsigned int charges = 0;
    double energies = 0.;
    if(!deposit_position_.empty()) {
        // Prepare charge deposits for this event
        std::vector<DepositedCharge> deposits;
@@ -224,6 +234,10 @@ void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* mess
            charges += 2 * charge;
            total_deposited_charge_ += 2 * charge;

            auto edep = deposit_energy_.at(i);
            total_deposited_energy_ += edep;
            energies += edep;

            // Match deposit with mc particle if possible
            auto track_id = deposit_to_id_.at(i);

@@ -251,10 +265,12 @@ void SensitiveDetectorActionG4::dispatchMessages(Module* module, Messenger* mess
    }
    // Store the number of charge carriers:
    deposited_charge_ = charges;
    deposited_energy_ = energies;

    // Clear deposit information for next event
    deposit_position_.clear();
    deposit_charge_.clear();
    deposit_energy_.clear();
    deposit_time_.clear();

    // Clear link tables for next event
+15 −0
Original line number Diff line number Diff line
@@ -62,6 +62,18 @@ namespace allpix {
         */
        unsigned int getDepositedCharge() const;

        /**
         * @brief Get total energy deposited in the sensitive device bound to this action
         */
        double getTotalDepositedEnergy() const;

        /**
         * @brief Get the energy deposited in the sensitive device for this event only.
         * @warning The correct number is only available after dispatching the message, before it refers to the previous
         * event.
         */
        double getDepositedEnergy() const;

        /**
         * @brief Get the name of the sensitive device bound to this action
         */
@@ -106,10 +118,13 @@ namespace allpix {
        // Statistics of total and per-event deposited charge
        unsigned int total_deposited_charge_{};
        unsigned int deposited_charge_{};
        double total_deposited_energy_{};
        double deposited_energy_{};

        // List of positions for deposits
        std::vector<ROOT::Math::XYZPoint> deposit_position_;
        std::vector<unsigned int> deposit_charge_;
        std::vector<double> deposit_energy_;
        std::vector<double> deposit_time_;

        // List of begin points for tracks