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

Rename Timepix3Clustering -> Clustering4D

parent 932197e6
......@@ -3,8 +3,7 @@ CORRYVRECKAN_DETECTOR_MODULE(MODULE_NAME)
# Add source files to library
CORRYVRECKAN_MODULE_SOURCES(${MODULE_NAME}
Timepix3Clustering.cpp
# ADD SOURCE FILES HERE...
Clustering4D.cpp
)
# Provide standard install target
......
#include "Timepix3Clustering.h"
#include "Clustering4D.h"
using namespace corryvreckan;
using namespace std;
Timepix3Clustering::Timepix3Clustering(Configuration config, std::shared_ptr<Detector> detector)
Clustering4D::Clustering4D(Configuration config, std::shared_ptr<Detector> detector)
: Module(std::move(config), detector), m_detector(detector) {
timingCut = m_config.get<double>("timingCut", static_cast<double>(Units::convert(100, "ns"))); // 100 ns
......@@ -11,7 +11,7 @@ Timepix3Clustering::Timepix3Clustering(Configuration config, std::shared_ptr<Det
neighbour_radius_col = m_config.get<int>("neighbour_radius_col", 1);
}
void Timepix3Clustering::initialise() {
void Clustering4D::initialise() {
// Cluster plots
std::string title = m_detector->name() + " Cluster size;cluster size;events";
......@@ -27,11 +27,11 @@ void Timepix3Clustering::initialise() {
}
// Sort function for pixels from low to high times
bool Timepix3Clustering::sortByTime(Pixel* pixel1, Pixel* pixel2) {
bool Clustering4D::sortByTime(Pixel* pixel1, Pixel* pixel2) {
return (pixel1->timestamp() < pixel2->timestamp());
}
StatusCode Timepix3Clustering::run(Clipboard* clipboard) {
StatusCode Clustering4D::run(Clipboard* clipboard) {
// Check if they are a Timepix3
if(m_detector->type() != "Timepix3") {
......@@ -132,7 +132,7 @@ StatusCode Timepix3Clustering::run(Clipboard* clipboard) {
}
// Check if a pixel touches any of the pixels in a cluster
bool Timepix3Clustering::touching(Pixel* neighbour, Cluster* cluster) {
bool Clustering4D::touching(Pixel* neighbour, Cluster* cluster) {
bool Touching = false;
......@@ -152,7 +152,7 @@ bool Timepix3Clustering::touching(Pixel* neighbour, Cluster* cluster) {
}
// Check if a pixel is close in time to the pixels of a cluster
bool Timepix3Clustering::closeInTime(Pixel* neighbour, Cluster* cluster) {
bool Clustering4D::closeInTime(Pixel* neighbour, Cluster* cluster) {
bool CloseInTime = false;
......@@ -166,7 +166,7 @@ bool Timepix3Clustering::closeInTime(Pixel* neighbour, Cluster* cluster) {
return CloseInTime;
}
void Timepix3Clustering::calculateClusterCentre(Cluster* cluster) {
void Clustering4D::calculateClusterCentre(Cluster* cluster) {
// Empty variables to calculate cluster position
double row(0), column(0), tot(0);
......
#ifndef TIMEPIX3CLUSTERING_H
#define TIMEPIX3CLUSTERING_H 1
#ifndef CLUSTERING4D_H
#define CLUSTERING4D_H 1
#include <TCanvas.h>
#include <TH1F.h>
......@@ -12,12 +12,12 @@
namespace corryvreckan {
/** @ingroup Modules
*/
class Timepix3Clustering : public Module {
class Clustering4D : public Module {
public:
// Constructors and destructors
Timepix3Clustering(Configuration config, std::shared_ptr<Detector> detector);
~Timepix3Clustering() {}
Clustering4D(Configuration config, std::shared_ptr<Detector> detector);
~Clustering4D() {}
// Functions
void initialise();
......@@ -42,4 +42,4 @@ namespace corryvreckan {
int neighbour_radius_col;
};
} // namespace corryvreckan
#endif // TIMEPIX3CLUSTERING_H
#endif // CLUSTERING4D_H
......@@ -39,7 +39,7 @@ OnlineMonitor::OnlineMonitor(Configuration config, std::vector<std::shared_ptr<D
canvas_cy2d =
m_config.getMatrix<std::string>("CorrelationY2D", {{"TestAlgorithm/%DETECTOR%/correlationY_2Dlocal", "colz"}});
canvas_charge = m_config.getMatrix<std::string>("ChargeDistributions", {{"Timepix3Clustering/%DETECTOR%/clusterTot"}});
canvas_charge = m_config.getMatrix<std::string>("ChargeDistributions", {{"Clustering4D/%DETECTOR%/clusterTot"}});
canvas_time = m_config.getMatrix<std::string>("EventTimes", {{"TestAlgorithm/%DETECTOR%/eventTimes"}});
}
......
#include "Timepix3Clustering.h"
using namespace corryvreckan;
using namespace std;
<<<<<<< HEAD
Timepix3Clustering::Timepix3Clustering(Configuration config, std::shared_ptr<Detector> detector)
: Module(std::move(config), std::move(detector)) {
timingCut = m_config.get<double>("timingCut", Units::convert(100, "ns")); // 100 ns
=======
Timepix3Clustering::Timepix3Clustering(Configuration config, std::vector<std::shared_ptr<Detector>> detectors)
: Module(std::move(config), std::move(detectors)) {
timingCut = m_config.get<double>("timingCut", static_cast<double>(Units::convert(100, "ns"))); // 100 ns
>>>>>>> compilerwarnings
neighbour_radius_row = m_config.get<int>("neighbour_radius_row", 1);
neighbour_radius_col = m_config.get<int>("neighbour_radius_col", 1);
}
void Timepix3Clustering::initialise() {
auto detector = get_detectors().front();
// Cluster plots
string name = "clusterSize_" + detector->name();
clusterSize = new TH1F(name.c_str(), name.c_str(), 100, 0, 100);
name = "clusterWidthRow_" + detector->name();
clusterWidthRow = new TH1F(name.c_str(), name.c_str(), 25, 0, 25);
name = "clusterWidthColumn_" + detector->name();
clusterWidthColumn = new TH1F(name.c_str(), name.c_str(), 100, 0, 100);
name = "clusterTot_" + detector->name();
clusterTot = new TH1F(name.c_str(), name.c_str(), 10000, 0, 100000);
name = "clusterPositionGlobal_" + detector->name();
clusterPositionGlobal = new TH2F(name.c_str(), name.c_str(), 400, -10., 10., 400, -10., 10.);
}
// Sort function for pixels from low to high times
bool Timepix3Clustering::sortByTime(Pixel* pixel1, Pixel* pixel2) {
return (pixel1->timestamp() < pixel2->timestamp());
}
StatusCode Timepix3Clustering::run(Clipboard* clipboard) {
auto detector = get_detectors().front();
// Check if they are a Timepix3
if(detector->type() != "Timepix3") {
return Success;
}
// Get the pixels
Pixels* pixels = (Pixels*)clipboard->get(detector->name(), "pixels");
if(pixels == NULL) {
LOG(DEBUG) << "Detector " << detector->name() << " does not have any pixels on the clipboard";
return Success;
}
LOG(DEBUG) << "Picked up " << pixels->size() << " pixels for device " << detector->name();
// Sort the pixels from low to high timestamp
std::sort(pixels->begin(), pixels->end(), sortByTime);
int totalPixels = pixels->size();
// Make the cluster storage
Clusters* deviceClusters = new Clusters();
// Keep track of which pixels are used
map<Pixel*, bool> used;
// Start to cluster
for(int iP = 0; iP < pixels->size(); iP++) {
Pixel* pixel = (*pixels)[iP];
// Check if pixel is used
if(used[pixel])
continue;
<<<<<<< HEAD
if(pixel->adc() == 0.)
continue;
// Make the new cluster object
Cluster* cluster = new Cluster();
LOG(DEBUG) << "==== New cluster";
// Keep adding hits to the cluster until no more are found
cluster->addPixel(pixel);
double clusterTime = pixel->timestamp();
used[pixel] = true;
LOG(DEBUG) << "Adding pixel: " << pixel->row() << "," << pixel->column();
int nPixels = 0;
while(cluster->size() != nPixels) {
nPixels = cluster->size();
// Loop over all pixels
for(int iNeighbour = (iP + 1); iNeighbour < totalPixels; iNeighbour++) {
Pixel* neighbour = (*pixels)[iNeighbour];
// Check if they are compatible in time with the cluster pixels
if((neighbour->timestamp() - clusterTime) > timingCut)
break;
// if(!closeInTime(neighbour,cluster)) break;
// Check if they have been used
if(used[neighbour])
continue;
if(neighbour->adc() == 0.)
continue;
// Check if they are touching cluster pixels
if(!touching(neighbour, cluster))
continue;
// Add to cluster
cluster->addPixel(neighbour);
clusterTime = neighbour->timestamp();
used[neighbour] = true;
LOG(DEBUG) << "Adding pixel: " << neighbour->row() << "," << neighbour->column() << " time "
<< Units::display(neighbour->timestamp(), {"ns", "us", "s"});
=======
// Get the pixels
Pixels* pixels = reinterpret_cast<Pixels*>(clipboard->get(detector->name(), "pixels"));
if(pixels == nullptr) {
LOG(DEBUG) << "Detector " << detector->name() << " does not have any pixels on the clipboard";
continue;
}
LOG(DEBUG) << "Picked up " << pixels->size() << " pixels for device " << detector->name();
// if(pixels->size() > 500.){
// LOG(DEBUG) <<"Skipping large event with "<<pixels->size()<<" pixels
// for device "<<detector->name();
// continue;
// }
// Sort the pixels from low to high timestamp
std::sort(pixels->begin(), pixels->end(), sortByTime);
size_t totalPixels = pixels->size();
// Make the cluster storage
Clusters* deviceClusters = new Clusters();
// Keep track of which pixels are used
map<Pixel*, bool> used;
// Start to cluster
for(size_t iP = 0; iP < pixels->size(); iP++) {
Pixel* pixel = (*pixels)[iP];
// Check if pixel is used
if(used[pixel]) {
continue;
}
if(pixel->adc() == 0.) {
continue;
}
// Make the new cluster object
Cluster* cluster = new Cluster();
LOG(DEBUG) << "==== New cluster";
// Keep adding hits to the cluster until no more are found
cluster->addPixel(pixel);
double clusterTime = pixel->timestamp();
used[pixel] = true;
LOG(DEBUG) << "Adding pixel: " << pixel->row() << "," << pixel->column();
size_t nPixels = 0;
while(cluster->size() != nPixels) {
nPixels = cluster->size();
// Loop over all pixels
for(size_t iNeighbour = (iP + 1); iNeighbour < totalPixels; iNeighbour++) {
Pixel* neighbour = (*pixels)[iNeighbour];
// Check if they are compatible in time with the cluster pixels
if((neighbour->timestamp() - clusterTime) > timingCut)
break;
// if(!closeInTime(neighbour,cluster)) break;
// Check if they have been used
if(used[neighbour])
continue;
if(neighbour->adc() == 0.)
continue;
// Check if they are touching cluster pixels
if(!touching(neighbour, cluster))
continue;
// Add to cluster
cluster->addPixel(neighbour);
clusterTime = neighbour->timestamp();
used[neighbour] = true;
LOG(DEBUG) << "Adding pixel: " << neighbour->row() << "," << neighbour->column() << " time "
<< Units::display(neighbour->timestamp(), {"ns", "us", "s"});
}
>>>>>>> compilerwarnings
}
}
// Finalise the cluster and save it
calculateClusterCentre(cluster);
<<<<<<< HEAD
// Fill cluster histograms
clusterSize->Fill(cluster->size());
clusterWidthRow->Fill(cluster->rowWidth());
clusterWidthColumn->Fill(cluster->columnWidth());
clusterTot->Fill(cluster->tot());
clusterPositionGlobal->Fill(cluster->globalX(), cluster->globalY());
=======
// Fill cluster histograms
clusterSize[detector->name()]->Fill(static_cast<double>(cluster->size()));
clusterWidthRow[detector->name()]->Fill(cluster->rowWidth());
clusterWidthColumn[detector->name()]->Fill(cluster->columnWidth());
clusterTot[detector->name()]->Fill(cluster->tot());
clusterPositionGlobal[detector->name()]->Fill(cluster->globalX(), cluster->globalY());
>>>>>>> compilerwarnings
deviceClusters->push_back(cluster);
}
<<<<<<< HEAD
// Put the clusters on the clipboard
if(deviceClusters->size() > 0) {
clipboard->put(detector->name(), "clusters", (Objects*)deviceClusters);
=======
// Put the clusters on the clipboard
if(deviceClusters->size() > 0) {
clipboard->put(detector->name(), "clusters", reinterpret_cast<Objects*>(deviceClusters));
}
LOG(DEBUG) << "Made " << deviceClusters->size() << " clusters for device " << detector->name();
>>>>>>> compilerwarnings
}
LOG(DEBUG) << "Made " << deviceClusters->size() << " clusters for device " << detector->name();
return Success;
}
// Check if a pixel touches any of the pixels in a cluster
bool Timepix3Clustering::touching(Pixel* neighbour, Cluster* cluster) {
bool Touching = false;
for(auto pixel : (*cluster->pixels())) {
int row_distance = abs(pixel->row() - neighbour->row());
int col_distance = abs(pixel->column() - neighbour->column());
if(row_distance <= neighbour_radius_row && col_distance <= neighbour_radius_col) {
if(row_distance > 1 || col_distance > 1) {
cluster->setSplit(true);
}
Touching = true;
break;
}
}
return Touching;
}
// Check if a pixel is close in time to the pixels of a cluster
bool Timepix3Clustering::closeInTime(Pixel* neighbour, Cluster* cluster) {
bool CloseInTime = false;
Pixels* pixels = cluster->pixels();
for(size_t iPix = 0; iPix < pixels->size(); iPix++) {
double timeDifference = abs(neighbour->timestamp() - (*pixels)[iPix]->timestamp());
if(timeDifference < timingCut)
CloseInTime = true;
}
return CloseInTime;
}
void Timepix3Clustering::calculateClusterCentre(Cluster* cluster) {
// Empty variables to calculate cluster position
double row(0), column(0), tot(0);
// Get the pixels on this cluster
Pixels* pixels = cluster->pixels();
string detectorID = (*pixels)[0]->detectorID();
double timestamp = (*pixels)[0]->timestamp();
// Loop over all pixels
for(auto& pixel : (*pixels)) {
double pixelToT = pixel->adc();
if(pixelToT == 0) {
LOG(DEBUG) << "Pixel with ToT 0!";
pixelToT = 1;
}
tot += pixelToT;
row += (pixel->row() * pixelToT);
column += (pixel->column() * pixelToT);
if(pixel->timestamp() < timestamp)
timestamp = pixel->timestamp();
}
// Row and column positions are tot-weighted
row /= (tot > 0 ? tot : 1);
column /= (tot > 0 ? tot : 1);
auto detector = get_detector(detectorID);
// Create object with local cluster position
PositionVector3D<Cartesian3D<double>> positionLocal(detector->pitch().X() * (column - detector->nPixelsX() / 2),
detector->pitch().Y() * (row - detector->nPixelsY() / 2),
0);
// Calculate global cluster position
PositionVector3D<Cartesian3D<double>> positionGlobal = detector->localToGlobal(positionLocal);
// Set the cluster parameters
cluster->setRow(row);
cluster->setColumn(column);
cluster->setTot(tot);
// Set uncertainty on position from intrinstic detector resolution:
cluster->setError(detector->resolution());
cluster->setTimestamp(timestamp);
cluster->setDetectorID(detectorID);
cluster->setClusterCentre(positionGlobal.X(), positionGlobal.Y(), positionGlobal.Z());
cluster->setClusterCentreLocal(positionLocal.X(), positionLocal.Y(), positionLocal.Z());
}
#ifndef TIMEPIX3CLUSTERING_H
#define TIMEPIX3CLUSTERING_H 1
#include <iostream>
#include <TCanvas.h>
#include <TH1F.h>
#include <TH2F.h>
#include "core/module/Module.hpp"
#include "objects/Cluster.h"
#include "objects/Pixel.h"
namespace corryvreckan {
/** @ingroup Modules
*/
class Timepix3Clustering : public Module {
public:
// Constructors and destructors
Timepix3Clustering(Configuration config, std::shared_ptr<Detector> detectors);
~Timepix3Clustering() {}
// Functions
void initialise();
StatusCode run(Clipboard* clipboard);
<<<<<<< HEAD
private:
=======
void finalise();
private:
static bool sortByTime(Pixel* pixel1, Pixel* pixel2);
>>>>>>> compilerwarnings
void calculateClusterCentre(Cluster*);
bool touching(Pixel*, Cluster*);
bool closeInTime(Pixel*, Cluster*);
// Cluster histograms
TH1F* clusterSize;
TH1F* clusterWidthRow;
TH1F* clusterWidthColumn;
TH1F* clusterTot;
TH2F* clusterPositionGlobal;
double timingCut;
int neighbour_radius_row;
int neighbour_radius_col;
};
} // namespace corryvreckan
#endif // TIMEPIX3CLUSTERING_H
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