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

Prealignment: make DETECTOR_MODULE and improve plots

parent 33c8774e
# Define module and return the generated name as MODULE_NAME
CORRYVRECKAN_GLOBAL_MODULE(MODULE_NAME)
CORRYVRECKAN_DETECTOR_MODULE(MODULE_NAME)
# Add source files to library
CORRYVRECKAN_MODULE_SOURCES(${MODULE_NAME}
......
......@@ -3,8 +3,9 @@
using namespace corryvreckan;
using namespace std;
Prealignment::Prealignment(Configuration config, std::vector<std::shared_ptr<Detector>> detectors)
: Module(std::move(config), std::move(detectors)) {
Prealignment::Prealignment(Configuration config, std::shared_ptr<Detector> detector)
: Module(std::move(config), detector), m_detector(detector) {
LOG(INFO) << "Starting prealignment of detectors";
max_correlation_rms = m_config.get<double>("max_correlation_rms", static_cast<double>(Units::convert(6.0, "mm")));
damping_factor = m_config.get<double>("damping_factor", 1.0);
......@@ -14,78 +15,71 @@ Prealignment::Prealignment(Configuration config, std::vector<std::shared_ptr<Det
}
void Prealignment::initialise() {
// Make histograms for each Timepix3
for(auto& detector : get_detectors()) {
LOG(INFO) << "Booking histograms for detector " << detector->name();
// get the reference detector:
std::shared_ptr<Detector> reference = get_reference();
// get the reference detector:
std::shared_ptr<Detector> reference = get_reference();
// Correlation plots
string name = "correlationX_" + detector->name();
correlationX[detector->name()] = new TH1F(name.c_str(), name.c_str(), 1000, -10., 10.);
name = "correlationY_" + detector->name();
correlationY[detector->name()] = new TH1F(name.c_str(), name.c_str(), 1000, -10., 10.);
// 2D correlation plots (pixel-by-pixel, local coordinates):
name = "correlationX_2Dlocal_" + detector->name();
correlationX2Dlocal[detector->name()] = new TH2F(name.c_str(),
name.c_str(),
detector->nPixelsX(),
0,
detector->nPixelsX(),
reference->nPixelsX(),
0,
reference->nPixelsX());
name = "correlationY_2Dlocal_" + detector->name();
correlationY2Dlocal[detector->name()] = new TH2F(name.c_str(),
name.c_str(),
detector->nPixelsY(),
0,
detector->nPixelsY(),
reference->nPixelsY(),
0,
reference->nPixelsY());
name = "correlationX_2D_" + detector->name();
correlationX2D[detector->name()] = new TH2F(name.c_str(), name.c_str(), 100, -10., 10., 100, -10., 10.);
name = "correlationY_2D_" + detector->name();
correlationY2D[detector->name()] = new TH2F(name.c_str(), name.c_str(), 100, -10., 10., 100, -10., 10.);
}
// Correlation plots
std::string 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.);
// 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(),
m_detector->nPixelsX(),
0,
m_detector->nPixelsX(),
reference->nPixelsX(),
0,
reference->nPixelsX());
title = m_detector->name() + ": 2D correlation Y (local);y [px];y_{ref} [px];events";
correlationY2Dlocal = new TH2F("correlationY_2Dlocal",
title.c_str(),
m_detector->nPixelsY(),
0,
m_detector->nPixelsY(),
reference->nPixelsY(),
0,
reference->nPixelsY());
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.);
}
StatusCode Prealignment::run(Clipboard* clipboard) {
// Loop over all Timepix3 and make plots
for(auto& detector : get_detectors()) {
// Get the clusters
Clusters* clusters = reinterpret_cast<Clusters*>(clipboard->get(detector->name(), "clusters"));
if(clusters == nullptr) {
LOG(DEBUG) << "Detector " << detector->name() << " does not have any clusters on the clipboard";
continue;
}
// 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";
return Success;
}
// 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";
continue;
}
// 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";
return Success;
}
// Loop over all clusters and fill histograms
for(auto& cluster : (*clusters)) {
// Loop over reference plane pixels to make correlation plots
for(auto& refCluster : (*referenceClusters)) {
double timeDifference = refCluster->timestamp() - cluster->timestamp();
// Loop over all clusters and fill histograms
for(auto& cluster : (*clusters)) {
// Loop over reference plane pixels to make correlation plots
for(auto& refCluster : (*referenceClusters)) {
double timeDifference = refCluster->timestamp() - cluster->timestamp();
// Correlation plots
if(abs(timeDifference) < timingCut) {
correlationX[detector->name()]->Fill(refCluster->globalX() - cluster->globalX());
correlationX2D[detector->name()]->Fill(cluster->globalX(), refCluster->globalX());
correlationX2Dlocal[detector->name()]->Fill(cluster->column(), refCluster->column());
correlationY[detector->name()]->Fill(refCluster->globalY() - cluster->globalY());
correlationY2D[detector->name()]->Fill(cluster->globalY(), refCluster->globalY());
correlationY2Dlocal[detector->name()]->Fill(cluster->row(), refCluster->row());
}
// Correlation plots
if(abs(timeDifference) < timingCut) {
correlationX->Fill(refCluster->globalX() - cluster->globalX());
correlationX2D->Fill(cluster->globalX(), refCluster->globalX());
correlationX2Dlocal->Fill(cluster->column(), refCluster->column());
correlationY->Fill(refCluster->globalY() - cluster->globalY());
correlationY2D->Fill(cluster->globalY(), refCluster->globalY());
correlationY2Dlocal->Fill(cluster->row(), refCluster->row());
}
}
}
......@@ -94,25 +88,25 @@ StatusCode Prealignment::run(Clipboard* clipboard) {
}
void Prealignment::finalise() {
for(auto& detector : get_detectors()) {
double rmsX = correlationX[detector->name()]->GetRMS();
double rmsY = correlationY[detector->name()]->GetRMS();
if(rmsX > max_correlation_rms or rmsY > max_correlation_rms) {
LOG(ERROR) << "Detector " << detector->name() << ": RMS is too wide for prealignment shifts";
LOG(ERROR) << "Detector " << detector->name() << ": RMS X = " << Units::display(rmsX, {"mm", "um"})
<< " , RMS Y = " << Units::display(rmsY, {"mm", "um"});
}
if(detector->role() != DetectorRole::REFERENCE) {
double mean_X = correlationX[detector->name()]->GetMean();
double mean_Y = correlationY[detector->name()]->GetMean();
LOG(INFO) << "Detector " << detector->name() << ": x = " << Units::display(mean_X, {"mm", "um"})
<< " , y = " << Units::display(mean_Y, {"mm", "um"});
LOG(INFO) << "Move in x by = " << Units::display(mean_X * damping_factor, {"mm", "um"})
<< " , and in y by = " << Units::display(mean_Y * damping_factor, {"mm", "um"});
double x = detector->displacement().X();
double y = detector->displacement().Y();
detector->displacementX(x + damping_factor * mean_X);
detector->displacementY(y + damping_factor * mean_Y);
}
double rmsX = correlationX->GetRMS();
double rmsY = correlationY->GetRMS();
if(rmsX > max_correlation_rms or rmsY > max_correlation_rms) {
LOG(ERROR) << "Detector " << m_detector->name() << ": RMS is too wide for prealignment shifts";
LOG(ERROR) << "Detector " << m_detector->name() << ": RMS X = " << Units::display(rmsX, {"mm", "um"})
<< " , RMS Y = " << Units::display(rmsY, {"mm", "um"});
}
if(m_detector->role() != DetectorRole::REFERENCE) {
double mean_X = correlationX->GetMean();
double mean_Y = correlationY->GetMean();
LOG(INFO) << "Detector " << m_detector->name() << ": x = " << Units::display(mean_X, {"mm", "um"})
<< " , y = " << Units::display(mean_Y, {"mm", "um"});
LOG(INFO) << "Move in x by = " << Units::display(mean_X * damping_factor, {"mm", "um"})
<< " , and in y by = " << Units::display(mean_Y * damping_factor, {"mm", "um"});
double x = m_detector->displacement().X();
double y = m_detector->displacement().Y();
m_detector->displacementX(x + damping_factor * mean_X);
m_detector->displacementY(y + damping_factor * mean_Y);
}
}
......@@ -16,7 +16,7 @@ namespace corryvreckan {
public:
// Constructors and destructors
Prealignment(Configuration config, std::vector<std::shared_ptr<Detector>> detectors);
Prealignment(Configuration config, std::shared_ptr<Detector> detector);
~Prealignment() {}
// Functions
......@@ -24,13 +24,16 @@ namespace corryvreckan {
StatusCode run(Clipboard* clipboard);
void finalise();
private:
std::shared_ptr<Detector> m_detector;
// Correlation plots
std::map<std::string, TH1F*> correlationX;
std::map<std::string, TH1F*> correlationY;
std::map<std::string, TH2F*> correlationX2Dlocal;
std::map<std::string, TH2F*> correlationY2Dlocal;
std::map<std::string, TH2F*> correlationX2D;
std::map<std::string, TH2F*> correlationY2D;
TH1F* correlationX;
TH1F* correlationY;
TH2F* correlationX2Dlocal;
TH2F* correlationY2Dlocal;
TH2F* correlationX2D;
TH2F* correlationY2D;
// Parameters which can be set by user
double max_correlation_rms;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment