Commit 8f1b1660 authored by Daniel Hynds's avatar Daniel Hynds
Browse files

use generic classes in timepix3 algorithms


Former-commit-id: 17f9c4f1bb9501c69d55fe2d90556f949e68414f
parent 5f2acdac
#include "BasicTracking.h"
#include "KDTreeTimepix3.h"
#include "KDTree.h"
#include "TCanvas.h"
BasicTracking::BasicTracking(bool debugging)
: Algorithm("BasicTracking"){
debug = debugging;
debug = false;
// Default values for cuts
timingCut = 200./1000000000.; // 200 ns
......@@ -37,18 +37,20 @@ void BasicTracking::initialise(Parameters* par){
residualsY[detectorID] = new TH1F(name.c_str(),name.c_str(),400,-0.2,0.2);
}
nTracksTotal = 0.;
}
StatusCode BasicTracking::run(Clipboard* clipboard){
if(debug) tcout<<"Start of event"<<endl;
// Container for all clusters, and detectors in tracking
// map<string,Timepix3Clusters*> clusters;
map<string,KDTreeTimepix3> trees;
map<string,KDTree> trees;
vector<string> detectors;
Timepix3Clusters* referenceClusters;
Clusters* referenceClusters;
// Output track container
Timepix3Tracks* tracks = new Timepix3Tracks();
Tracks* tracks = new Tracks();
// Loop over all Timepix3 and get clusters
for(int det = 0; det<parameters->nDetectors; det++){
......@@ -58,14 +60,16 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
if(parameters->detector[detectorID]->type() != "Timepix3") continue;
// Get the clusters
Timepix3Clusters* tempClusters = (Timepix3Clusters*)clipboard->get(detectorID,"clusters");
Clusters* tempClusters = (Clusters*)clipboard->get(detectorID,"clusters");
if(tempClusters == NULL){
if(debug) tcout<<"Detector "<<detectorID<<" does not have any clusters on the clipboard"<<endl;
}else{
// Store them
// clusters[detectorID] = tempClusters;
if(debug) tcout<<"Picked up "<<tempClusters->size()<<" clusters from "<<detectorID<<endl;
if(detectorID == parameters->reference) referenceClusters = tempClusters;
KDTreeTimepix3 clusterTree(*tempClusters);
KDTree clusterTree;
clusterTree.buildTimeTree(*tempClusters);
trees[detectorID] = clusterTree;
detectors.push_back(detectorID);
}
......@@ -76,7 +80,7 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
// Use the first plane as a seeding plane. For something quick, look a cluster in < 200 ns in the next plane, and continue
string reference = parameters->reference;
map<Timepix3Cluster*, bool> used;
map<Cluster*, bool> used;
// If no clusters on reference plane, stop
if(trees.count(reference) == 0) return Success;
......@@ -86,9 +90,10 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
for(int iSeedCluster=0;iSeedCluster<nSeedClusters;iSeedCluster++){
// Make a new track
Timepix3Track* track = new Timepix3Track();
if(debug) tcout<<"Looking at seed cluster "<<iSeedCluster<<endl;
Track* track = new Track();
// Get the cluster
Timepix3Cluster* cluster = (*referenceClusters)[iSeedCluster];
Cluster* cluster = (*referenceClusters)[iSeedCluster];
// Add the cluster to the track
track->addCluster(cluster);
track->setTimestamp(cluster->timestamp());
......@@ -105,7 +110,7 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
// // If the detector is excluded from tracking ignore it
// if(parameters->excludedFromTracking.count(detectorID) != 0) continue;
// // Get the closest cluster
// Timepix3Cluster* newCluster = getNearestCluster(cluster, used, clusters[detectorID]);
// Cluster* newCluster = getNearestCluster(cluster, used, clusters[detectorID]);
// if(newCluster == NULL) continue;
// // Add the cluster to the track
// track->addCluster(newCluster);
......@@ -119,7 +124,7 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
if(clusters[detectors[det]] == NULL) continue;
// If excluded from tracking ignore this plane
if(parameters->excludedFromTracking.count(detectors[det]) != 0) continue;
Timepix3Cluster* newCluster = getNearestCluster(timestamp, (*clusters[detectors[det]]) );
Cluster* newCluster = getNearestCluster(timestamp, (*clusters[detectors[det]]) );
if( ((newCluster->timestamp() - timestamp) / (4096.*40000000.)) > (10./1000000000.) ) continue;
// Check if spatially more than 200 um
if( abs(cluster->globalX() - newCluster->globalX()) > spatialCut || abs(cluster->globalY() - newCluster->globalY()) > spatialCut ) continue;
......@@ -135,12 +140,15 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
if(parameters->excludedFromTracking.count(detectors[det]) != 0) continue;
// Get all neighbours within 200 ns
Timepix3Cluster* closestCluster = NULL; double closestClusterDistance = spatialCut;
if(debug) tcout<<"Searching for neighbouring cluster on "<<detectors[det]<<endl;
if(debug) tcout<<"- cluster time is "<<cluster->timestamp()<<endl;
Cluster* closestCluster = NULL; double closestClusterDistance = spatialCut;
// tcout<<"About to get nearest neighbours in window: "<<timingCut<<endl;
Timepix3Clusters neighbours = trees[detectors[det]].getAllClustersInTimeWindow(cluster,timingCut);
Clusters neighbours = trees[detectors[det]].getAllClustersInTimeWindow(cluster,timingCut);
if(debug) tcout<<"- found "<<neighbours.size()<<" neighbours"<<endl;
for(int ne=0;ne<neighbours.size();ne++){
Timepix3Cluster* newCluster = neighbours[ne];
Cluster* newCluster = neighbours[ne];
double distance = sqrt((cluster->globalX() - newCluster->globalX())*(cluster->globalX() - newCluster->globalX()) + (cluster->globalY() - newCluster->globalY())*(cluster->globalY() - newCluster->globalY()));
......@@ -154,6 +162,7 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
if(closestCluster == NULL) continue;
// Add the cluster to the track
if(debug) tcout<<"- added cluster to track"<<endl;
track->addCluster(closestCluster);
}//*/
......@@ -176,9 +185,9 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
trackAngleY->Fill(atan(track->m_direction.Y()));
// Make residuals
Timepix3Clusters trackClusters = track->clusters();
Clusters trackClusters = track->clusters();
for(int iTrackCluster=0; iTrackCluster<trackClusters.size(); iTrackCluster++){
Timepix3Cluster* trackCluster = trackClusters[iTrackCluster];
Cluster* trackCluster = trackClusters[iTrackCluster];
string detectorID = trackCluster->detectorID();
ROOT::Math::XYZPoint intercept = track->intercept(trackCluster->globalZ());
residualsX[detectorID]->Fill(intercept.X() - trackCluster->globalX());
......@@ -189,25 +198,29 @@ StatusCode BasicTracking::run(Clipboard* clipboard){
// Save the tracks on the clipboard
if(tracks->size() > 0){
clipboard->put("Timepix3","tracks",(TestBeamObjects*)tracks);
clipboard->put("tracks",(TestBeamObjects*)tracks);
tracksPerEvent->Fill(tracks->size());
}
nTracksTotal+=tracks->size();
// cout<<", produced "<<nTracksTotal<<" tracks";
// Clean up tree objects
// for(int det = 0; det<parameters->nDetectors; det++){
// string detectorID = parameters->detectors[det];
// if(trees.count(detectorID) != 0) delete trees[detectorID];
// }
if(debug) tcout<<"End of event"<<endl;
return Success;
}
Timepix3Cluster* BasicTracking::getNearestCluster(long long int timestamp, Timepix3Clusters clusters){
Cluster* BasicTracking::getNearestCluster(long long int timestamp, Clusters clusters){
Timepix3Cluster* bestCluster = NULL;
Cluster* bestCluster = NULL;
// Loop over all clusters and return the one with the closest timestamp
for(int iCluster=0;iCluster<clusters.size();iCluster++){
Timepix3Cluster* cluster = clusters[iCluster];
Cluster* cluster = clusters[iCluster];
if(bestCluster == NULL) bestCluster = cluster;
if(abs(cluster->timestamp() - timestamp) < abs(bestCluster->timestamp()-timestamp)) bestCluster = cluster;
}
......
......@@ -6,9 +6,9 @@
#include "TH1F.h"
#include "TH2F.h"
#include "TCanvas.h"
#include "Timepix3Pixel.h"
#include "Timepix3Cluster.h"
#include "Timepix3Track.h"
#include "Pixel.h"
#include "Cluster.h"
#include "Track.h"
class BasicTracking : public Algorithm {
......@@ -22,8 +22,8 @@ public:
StatusCode run(Clipboard*);
void finalise();
// Timepix3Cluster* getNearestCluster(Timepix3Cluster*, map<Timepix3Cluster*, bool>, Timepix3Clusters*);
Timepix3Cluster* getNearestCluster(long long int, Timepix3Clusters);
// Cluster* getNearestCluster(Cluster*, map<Cluster*, bool>, Clusters*);
Cluster* getNearestCluster(long long int, Clusters);
// Member variables
......@@ -41,6 +41,7 @@ public:
double timingCut;
double spatialCut;
int minHitsOnTrack;
double nTracksTotal;
};
......
......@@ -15,7 +15,7 @@ void Timepix3Clustering::initialise(Parameters* par){
}
// Sort function for pixels from low to high times
bool sortByTime(Timepix3Pixel* pixel1, Timepix3Pixel* pixel2){
bool sortByTime(Pixel* pixel1, Pixel* pixel2){
return (pixel1->m_timestamp < pixel2->m_timestamp);
}
......@@ -29,7 +29,7 @@ StatusCode Timepix3Clustering::run(Clipboard* clipboard){
if(parameters->detector[detectorID]->type() != "Timepix3") continue;
// Get the pixels
Timepix3Pixels* pixels = (Timepix3Pixels*)clipboard->get(detectorID,"pixels");
Pixels* pixels = (Pixels*)clipboard->get(detectorID,"pixels");
if(pixels == NULL){
if(debug) tcout<<"Detector "<<detectorID<<" does not have any pixels on the clipboard"<<endl;
continue;
......@@ -46,20 +46,20 @@ StatusCode Timepix3Clustering::run(Clipboard* clipboard){
int totalPixels = pixels->size();
// Make the cluster storage
Timepix3Clusters* deviceClusters = new Timepix3Clusters();
Clusters* deviceClusters = new Clusters();
// Keep track of which pixels are used
map<Timepix3Pixel*, bool> used;
map<Pixel*, bool> used;
// Start to cluster
for(int iP=0;iP<pixels->size();iP++){
Timepix3Pixel* pixel = (*pixels)[iP];
Pixel* pixel = (*pixels)[iP];
// Check if pixel is used
if(used[pixel]) continue;
// Make the new cluster object
Timepix3Cluster* cluster = new Timepix3Cluster();
Cluster* cluster = new Cluster();
if(debug) tcout<<"==== New cluster"<<endl;
// Keep adding hits to the cluster until no more are found
......@@ -73,7 +73,7 @@ StatusCode Timepix3Clustering::run(Clipboard* clipboard){
nPixels = cluster->size();
// Loop over all pixels
for(int iNeighbour=(iP+1);iNeighbour<totalPixels;iNeighbour++){
Timepix3Pixel* neighbour = (*pixels)[iNeighbour];
Pixel* neighbour = (*pixels)[iNeighbour];
// Check if they are compatible in time with the cluster pixels
if( (neighbour->m_timestamp - clusterTime) > timingCutInt ) break;
// if(!closeInTime(neighbour,cluster)) break;
......@@ -106,10 +106,10 @@ StatusCode Timepix3Clustering::run(Clipboard* clipboard){
}
// Check if a pixel touches any of the pixels in a cluster
bool Timepix3Clustering::touching(Timepix3Pixel* neighbour,Timepix3Cluster* cluster){
bool Timepix3Clustering::touching(Pixel* neighbour,Cluster* cluster){
bool Touching = false;
Timepix3Pixels pixels = cluster->pixels();
Pixels pixels = cluster->pixels();
for(int iPix=0;iPix<pixels.size();iPix++){
if( abs(pixels[iPix]->m_row - neighbour->m_row) <= 1 &&
......@@ -122,11 +122,11 @@ bool Timepix3Clustering::touching(Timepix3Pixel* neighbour,Timepix3Cluster* clus
}
// Check if a pixel is close in time to the pixels of a cluster
bool Timepix3Clustering::closeInTime(Timepix3Pixel* neighbour,Timepix3Cluster* cluster){
bool Timepix3Clustering::closeInTime(Pixel* neighbour,Cluster* cluster){
bool CloseInTime = false;
Timepix3Pixels pixels = cluster->pixels();
Pixels pixels = cluster->pixels();
for(int iPix=0;iPix<pixels.size();iPix++){
long long int timeDifference = abs(neighbour->m_timestamp - pixels[iPix]->m_timestamp);
......@@ -137,14 +137,14 @@ bool Timepix3Clustering::closeInTime(Timepix3Pixel* neighbour,Timepix3Cluster* c
}
void Timepix3Clustering::calculateClusterCentre(Timepix3Cluster* cluster){
void Timepix3Clustering::calculateClusterCentre(Cluster* cluster){
// Empty variables to calculate cluster position
double row(0), column(0), tot(0);
long long int timestamp;
// Get the pixels on this cluster
Timepix3Pixels pixels = cluster->pixels();
Pixels pixels = cluster->pixels();
string detectorID = pixels[0]->m_detectorID;
timestamp = pixels[0]->m_timestamp;
......
#ifndef TIMEPIX3CLUSTERING_H
#define TIMEPIX3CLUSTERING_H 1
#include "Timepix3Pixel.h"
#include "Timepix3Cluster.h"
#include "Pixel.h"
#include "Cluster.h"
#include "Algorithm.h"
#include <iostream>
#include "TH1F.h"
......@@ -20,9 +20,9 @@ public:
void initialise(Parameters*);
StatusCode run(Clipboard*);
void finalise();
void calculateClusterCentre(Timepix3Cluster*);
bool touching(Timepix3Pixel*,Timepix3Cluster*);
bool closeInTime(Timepix3Pixel*,Timepix3Cluster*);
void calculateClusterCentre(Cluster*);
bool touching(Pixel*,Cluster*);
bool closeInTime(Pixel*,Cluster*);
double timingCut;
long long int timingCutInt;
......
......@@ -77,7 +77,7 @@ void Timepix3EventLoader::initialise(Parameters* par){
tcout<<"Registering detector "<<detectorID<<endl;
if(parameters->detector.count(detectorID) == 0){tcout<<"Detector "<<detectorID<<" has no alignment/conditions loaded. Please check that it is in the alignment file"<<endl; return;}
parameters->registerDetector(detectorID);
// parameters->registerDetector(detectorID);
// Now that we have all of the data files and mask files for this detector, pass the mask file to parameters
tcout<<"Set mask file "<<trimdacfile<<endl;
......@@ -98,7 +98,8 @@ StatusCode Timepix3EventLoader::run(Clipboard* clipboard){
// loading a fixed number of pixels (ie. 2000 at a time)
int endOfFiles = 0; int devices = 0; int loadedData = 0;
cout<<"\rCurrent time: "<<parameters->currentTime<<flush;
// Loop through all registered detectors
for(int det = 0; det<parameters->nDetectors; det++){
......@@ -107,7 +108,7 @@ StatusCode Timepix3EventLoader::run(Clipboard* clipboard){
if(parameters->detector[detectorID]->type() != "Timepix3") continue;
// Make a new container for the data
Timepix3Pixels* deviceData = new Timepix3Pixels();
Pixels* deviceData = new Pixels();
SpidrSignals* spidrData = new SpidrSignals();
// Load the next chunk of data
......@@ -166,7 +167,7 @@ void Timepix3EventLoader::maskPixels(string detectorID, string trimdacfile){
}
// Function to load data for a given device, into the relevant container
bool Timepix3EventLoader::loadData(string detectorID, Timepix3Pixels* devicedata, SpidrSignals* spidrData){
bool Timepix3EventLoader::loadData(string detectorID, Pixels* devicedata, SpidrSignals* spidrData){
// if(detectorID == "W0019_F07") debug = true;
// if(detectorID != "W0019_F07") debug = false;
......@@ -337,7 +338,7 @@ bool Timepix3EventLoader::loadData(string detectorID, Timepix3Pixels* devicedata
}
// Otherwise create a new pixel object
Timepix3Pixel* pixel = new Timepix3Pixel(detectorID,row,col,(int)tot,time);
Pixel* pixel = new Pixel(detectorID,row,col,(int)tot,time);
devicedata->push_back(pixel);
npixels++;
......
......@@ -2,7 +2,7 @@
#define TIMEPIX3EVENTLOADER_H 1
#include "Algorithm.h"
#include "Timepix3Pixel.h"
#include "Pixel.h"
#include "SpidrSignal.h"
#include <stdio.h>
......@@ -19,7 +19,7 @@ public:
void finalise();
// Extra functions
bool loadData(string, Timepix3Pixels*, SpidrSignals*);
bool loadData(string, Pixels*, SpidrSignals*);
void maskPixels(string, string);
// Member variables
......
......@@ -39,7 +39,7 @@ StatusCode Timepix3MaskCreator::run(Clipboard* clipboard){
if(parameters->detector[detectorID]->type() != "Timepix3") continue;
// Get the pixels
Timepix3Pixels* pixels = (Timepix3Pixels*)clipboard->get(detectorID,"pixels");
Pixels* pixels = (Pixels*)clipboard->get(detectorID,"pixels");
if(pixels == NULL){
tcout<<"Detector "<<detectorID<<" does not have any pixels on the clipboard"<<endl;
return Failure;
......@@ -48,7 +48,7 @@ StatusCode Timepix3MaskCreator::run(Clipboard* clipboard){
// Loop over all pixels
for(int iP=0;iP<pixels->size();iP++){
Timepix3Pixel* pixel = (*pixels)[iP];
Pixel* pixel = (*pixels)[iP];
// Enter another pixel hit for this channel
int channelID = pixel->m_row + 256*pixel->m_column;
......
#ifndef TIMEPIX3MASKCREATOR_H
#define TIMEPIX3MASKCREATOR_H 1
#include "Timepix3Pixel.h"
#include "Pixel.h"
#include "Algorithm.h"
#include <iostream>
#include "TH1F.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