TestAlgorithm.cpp 6.51 KB
Newer Older
1
2
#include "TestAlgorithm.h"

3
using namespace corryvreckan;
4
using namespace std;
5

6
7
8
TestAlgorithm::TestAlgorithm(Configuration config, std::shared_ptr<Detector> detector)
    : Module(std::move(config), detector), m_detector(detector) {

9
10
11
    makeCorrelations = m_config.get<bool>("make_correlations", false);
    timingCut = m_config.get<double>("timing_cut", static_cast<double>(Units::convert(100, "ns")));
    do_timing_cut_ = m_config.get<bool>("do_timing_cut", false);
12
13
}

14
void TestAlgorithm::initialise() {
15
16
17
18
19
20
    LOG(DEBUG) << "Booking histograms for detector " << m_detector->name();

    // get the reference detector:
    std::shared_ptr<Detector> reference = get_reference();

    // Simple hit map
21
    std::string title = m_detector->name() + ": hitmap;x [px];y [px];events";
22
    hitmap = new TH2F("hitmap",
23
                      title.c_str(),
24
                      m_detector->nPixels().X(),
25
                      0,
26
27
                      m_detector->nPixels().X(),
                      m_detector->nPixels().Y(),
28
                      0,
29
                      m_detector->nPixels().Y());
30

31
32
33
34
35
36
    if(makeCorrelations) {
        // Correlation plots
        title = m_detector->name() + ": correlation X;x_{ref}-x [mm];events";
        correlationX = new TH1F("correlationX", title.c_str(), 1000, -10., 10.);
        title = m_detector->name() + ": correlation Y;y_{ref}-y [mm];events";
        correlationY = new TH1F("correlationY", title.c_str(), 1000, -10., 10.);
37

38
39
        // time correlation plot range should cover length of events. nanosecond binning.
        title = m_detector->name() + "Reference cluster time stamp - cluster time stamp;t_{ref}-t [ns];events";
40
41
        correlationTime =
            new TH1F("correlationTime", title.c_str(), static_cast<int>(2. * timingCut), -1 * timingCut, timingCut);
42
43
        title = m_detector->name() + "Reference cluster time stamp - cluster time stamp;t_{ref}-t [1/40MHz];events";
        correlationTimeInt = new TH1F("correlationTimeInt", title.c_str(), 8000, -40000, 40000);
44

45
46
47
48
        // 2D correlation plots (pixel-by-pixel, local coordinates):
        title = m_detector->name() + ": 2D correlation X (local);x [px];x_{ref} [px];events";
        correlationX2Dlocal = new TH2F("correlationX_2Dlocal",
                                       title.c_str(),
49
                                       m_detector->nPixels().X(),
50
                                       0,
51
52
                                       m_detector->nPixels().X(),
                                       reference->nPixels().X(),
53
                                       0,
54
                                       reference->nPixels().X());
55
56
57
        title = m_detector->name() + ": 2D correlation Y (local);y [px];y_{ref} [px];events";
        correlationY2Dlocal = new TH2F("correlationY_2Dlocal",
                                       title.c_str(),
58
                                       m_detector->nPixels().Y(),
59
                                       0,
60
61
                                       m_detector->nPixels().Y(),
                                       reference->nPixels().Y(),
62
                                       0,
63
                                       reference->nPixels().Y());
64

65
66
67
68
69
        title = m_detector->name() + ": 2D correlation X (global);x [mm];x_{ref} [mm];events";
        correlationX2D = new TH2F("correlationX_2D", title.c_str(), 100, -10., 10., 100, -10., 10.);
        title = m_detector->name() + ": 2D correlation Y (global);y [mm];y_{ref} [mm];events";
        correlationY2D = new TH2F("correlationY_2D", title.c_str(), 100, -10., 10., 100, -10., 10.);
    }
70
    // Timing plots
71
72
    title = m_detector->name() + ": event time;t [s];events";
    eventTimes = new TH1F("eventTimes", title.c_str(), 3000000, 0, 300);
73
74
}

75
StatusCode TestAlgorithm::run(std::shared_ptr<Clipboard> clipboard) {
Simon Spannagel's avatar
Simon Spannagel committed
76

77
78
79
80
    // Get the pixels
    Pixels* pixels = reinterpret_cast<Pixels*>(clipboard->get(m_detector->name(), "pixels"));
    if(pixels == nullptr) {
        LOG(DEBUG) << "Detector " << m_detector->name() << " does not have any pixels on the clipboard";
81
        return StatusCode::Success;
82
    }
Simon Spannagel's avatar
Simon Spannagel committed
83

84
85
86
87
88
89
90
    // Loop over all pixels and make hitmaps
    for(auto& pixel : (*pixels)) {
        // Hitmap
        hitmap->Fill(pixel->column(), pixel->row());
        // Timing plots
        eventTimes->Fill(static_cast<double>(Units::convert(pixel->timestamp(), "s")));
    }
Simon Spannagel's avatar
Simon Spannagel committed
91

92
93
94
95
    // Get the clusters
    Clusters* clusters = reinterpret_cast<Clusters*>(clipboard->get(m_detector->name(), "clusters"));
    if(clusters == nullptr) {
        LOG(DEBUG) << "Detector " << m_detector->name() << " does not have any clusters on the clipboard";
96
        return StatusCode::Success;
97
    }
Simon Spannagel's avatar
Simon Spannagel committed
98

99
100
101
102
103
    // Get clusters from reference detector
    auto reference = get_reference();
    Clusters* referenceClusters = reinterpret_cast<Clusters*>(clipboard->get(reference->name(), "clusters"));
    if(referenceClusters == nullptr) {
        LOG(DEBUG) << "Reference detector " << reference->name() << " does not have any clusters on the clipboard";
104
        return StatusCode::Success;
105
    }
Simon Spannagel's avatar
Simon Spannagel committed
106

107
108
109
110
111
112
113
114
115
116
117
    // Loop over all clusters and fill histograms
    if(makeCorrelations) {
        for(auto& cluster : (*clusters)) {
            // Loop over reference plane pixels to make correlation plots
            for(auto& refCluster : (*referenceClusters)) {

                double timeDifference = refCluster->timestamp() - cluster->timestamp();
                // in 40 MHz:
                long long int timeDifferenceInt = static_cast<long long int>(timeDifference / 25);

                // Correlation plots
118
                if(abs(timeDifference) < timingCut || !do_timing_cut_) {
119
120
                    correlationX->Fill(refCluster->global().x() - cluster->global().x());
                    correlationX2D->Fill(cluster->global().x(), refCluster->global().x());
121
                    correlationX2Dlocal->Fill(cluster->column(), refCluster->column());
122
                }
123
                if(abs(timeDifference) < timingCut || !do_timing_cut_) {
124
125
                    correlationY->Fill(refCluster->global().y() - cluster->global().y());
                    correlationY2D->Fill(cluster->global().y(), refCluster->global().y());
126
127
128
129
130
                    correlationY2Dlocal->Fill(cluster->row(), refCluster->row());
                }
                //                    correlationTime[m_detector->name()]->Fill(Units::convert(timeDifference, "s"));
                correlationTime->Fill(timeDifference); // time difference in ns
                correlationTimeInt->Fill(static_cast<double>(timeDifferenceInt));
131
            }
Simon Spannagel's avatar
Simon Spannagel committed
132
133
        }
    }
134

135
    return StatusCode::Success;
136
}