Commit f217b6df authored by Simon Spannagel's avatar Simon Spannagel
Browse files

Merge branch 'energy_deposit_histo' into 'master'

DepositionG4: add histogram diplaying the deposited energy per event

See merge request allpix-squared/allpix-squared!700
parents f98a871a 0baad526
Loading
Loading
Loading
Loading
+29 −6
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
@@ -419,17 +425,34 @@ void DepositionGeant4Module::construct_sensitive_detectors_and_fields(double fan
            LOG(TRACE) << "Creating output plots for detector " << sensitive_detector_action->getName();

            // Plot axis are in kilo electrons - convert from framework units!
            int maximum = static_cast<int>(Units::convert(config_.get<int>("output_plots_scale"), "ke"));
            int nbins = 5 * maximum;
            int maximum_charge = static_cast<int>(Units::convert(config_.get<int>("output_plots_scale"), "ke"));
            double maximum_energy =
                (static_cast<int>(maximum_charge / 2. * Units::convert(charge_creation_energy, "eV")) / 10) * 10 + 10;
            int nbins = 5 * maximum_charge;

            // 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);
                    charge_per_event_[sensitive_detector_action->getName()] =
                        CreateHistogram<TH1D>(plot_name.c_str(),
                                              "deposited charge per event;deposited charge [ke];events",
                                              nbins,
                                              0,
                                              maximum_charge);
                }

                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",
                                              nbins,
                                              0,
                                              maximum_energy);
                }
            }
        }
+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