fixes with beam data. End of file management, checks for no data, clearing of...

fixes with beam data. End of file management, checks for no data, clearing of objects from clipboard


Former-commit-id: aa5dad9cd50587cc46578062c65a2d7578c191ae

branches/trunk/algorithms/BasicTracking.C

@@ -80,15 +80,16 @@ int BasicTracking::run(Clipboard* clipboard){
     // Loop over each subsequent plane and look for a cluster within 100 ns
     for(int det=0; det<detectors.size(); det++){
       if(detectors[det] == reference) continue;
+      if(clusters[detectors[det]] == NULL) continue;
       Timepix3Cluster* 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()) > 0.3 || abs(cluster->globalY() - newCluster->globalY()) > 0.3 ) continue;
+      if( abs(cluster->globalX() - newCluster->globalX()) > 0.2 || abs(cluster->globalY() - newCluster->globalY()) > 0.2 ) continue;
       track->addCluster(newCluster);
     }
    
     // Now should have a track with one cluster from each plane
-    if(track->nClusters() < 5){
+    if(track->nClusters() < 6){
       delete track;
       continue;
     }

branches/trunk/algorithms/GUI.C

@@ -22,6 +22,8 @@ void startDisplay(void* gui){
 //  ((GUI*)gui)->trackCanvas->Divide();
   ((GUI*)gui)->hitmapCanvas = new TCanvas("HitMapCanvas","Hit map canvas");
   ((GUI*)gui)->hitmapCanvas->Divide(ceil(nDetectors/2.),2);
+  ((GUI*)gui)->globalHitmapCanvas = new TCanvas("GlobalHitmapCanvas","Global hit map canvas");
+  ((GUI*)gui)->globalHitmapCanvas->Divide(ceil(nDetectors/2.),2);
   ((GUI*)gui)->residualsCanvas = new TCanvas("ResidualsCanvas","Residuals canvas");
   ((GUI*)gui)->residualsCanvas->Divide(ceil(nDetectors/2.),2);
 
@@ -50,6 +52,8 @@ void GUI::initialise(Parameters* par){
     hitmap[detectorID] = (TH2F*)gDirectory->Get(hitmapHisto.c_str());
     string residualHisto = "/tbAnalysis/BasicTracking/residualsX_"+detectorID;
     residuals[detectorID] = (TH1F*)gDirectory->Get(residualHisto.c_str());
+    string globalHitmapHisto = "/tbAnalysis/TestAlgorithm/clusterPositionGlobal_"+detectorID;
+    globalHitmap[detectorID] = (TH2F*)gDirectory->Get(globalHitmapHisto.c_str());
   }
 
   // Set event counter
@@ -92,7 +96,26 @@ int GUI::run(Clipboard* clipboard){
     hitmapCanvas->Paint();
     hitmapCanvas->Update();
   }
+
+  //-----------------------------------------
+  // Draw the objects on the globalHitmap canvas
+  //-----------------------------------------
   
+  // Update the canvas
+  if(eventNumber == 0){
+    for(int det = 0; det<parameters->nDetectors; det++){
+      string detectorID = parameters->detectors[det];
+      globalHitmap[detectorID]->SetTitle(detectorID.c_str());
+      globalHitmapCanvas->cd(det+1);
+      globalHitmap[detectorID]->Draw("colz");
+    }
+  }
+  // Update the canvas
+  if(eventNumber%updateNumber == 0) {
+    globalHitmapCanvas->Paint();
+    globalHitmapCanvas->Update();
+  }
+
   //-----------------------------------------
   // Draw the objects on the residuals canvas
   //-----------------------------------------

branches/trunk/algorithms/GUI.h

@@ -30,10 +30,12 @@ public:
   // Canvases to display plots
   TCanvas* trackCanvas;
   TCanvas* hitmapCanvas;
+  TCanvas* globalHitmapCanvas;
   TCanvas* residualsCanvas;
   
   // Plot holders
   map<string, TH2F*> hitmap;
+  map<string, TH2F*> globalHitmap;
   map<string, TH1F*> residuals;
   
   // Thread to allow display to run in separate thread

branches/trunk/algorithms/TestAlgorithm.C

@@ -49,13 +49,6 @@ void TestAlgorithm::initialise(Parameters* par){
 
 int TestAlgorithm::run(Clipboard* clipboard){
   
-  // Get clusters from reference detector
-  Timepix3Clusters* referenceClusters = (Timepix3Clusters*)clipboard->get(parameters->reference,"clusters");
-  if(referenceClusters == NULL){
-    tcout<<"Reference detector "<<parameters->reference<<" does not have any clusters on the clipboard"<<endl;
-    return 1;
-  }
-  
   // Loop over all Timepix3 and make plots
   for(int det = 0; det<parameters->nDetectors; det++){
     
@@ -67,14 +60,14 @@ int TestAlgorithm::run(Clipboard* clipboard){
     Timepix3Pixels* pixels = (Timepix3Pixels*)clipboard->get(detectorID,"pixels");
     if(pixels == NULL){
       tcout<<"Detector "<<detectorID<<" does not have any pixels on the clipboard"<<endl;
-      return 1;
+      continue;
     }
     
     // Get the clusters
     Timepix3Clusters* clusters = (Timepix3Clusters*)clipboard->get(detectorID,"clusters");
     if(clusters == NULL){
       tcout<<"Detector "<<detectorID<<" does not have any clusters on the clipboard"<<endl;
-      return 1;
+      continue;
     }
     
     // Loop over all pixels and make hitmaps
@@ -84,13 +77,20 @@ int TestAlgorithm::run(Clipboard* clipboard){
       Timepix3Pixel* pixel = (*pixels)[iP];
       
       // Hitmap
-      hitmap[detectorID]->Fill(pixel->m_row,pixel->m_column);
+      hitmap[detectorID]->Fill(pixel->m_column,pixel->m_row);
       
       // Timing plots
       eventTimes[detectorID]->Fill((double)pixel->m_timestamp / (4096.*40000000.) );
       
     }
     
+    // Get clusters from reference detector
+    Timepix3Clusters* referenceClusters = (Timepix3Clusters*)clipboard->get(parameters->reference,"clusters");
+    if(referenceClusters == NULL){
+      tcout<<"Reference detector "<<parameters->reference<<" does not have any clusters on the clipboard"<<endl;
+      continue;
+    }
+
     // Loop over all clusters and fill histograms
     for(int iCluster=0;iCluster<clusters->size();iCluster++){
 
@@ -111,8 +111,8 @@ int TestAlgorithm::run(Clipboard* clipboard){
         double timeDifference = (double)(refCluster->timestamp() - cluster->timestamp()) / (4096.*40000000.);
         
         // Correlation plots
-        if( abs(timeDifference) < 0.001 ) correlationX[detectorID]->Fill(refCluster->globalX() - cluster->globalX());
-        if( abs(timeDifference) < 0.001 ) correlationY[detectorID]->Fill(refCluster->globalY() - cluster->globalY());
+        if( abs(timeDifference) < 0.0001 ) correlationX[detectorID]->Fill(refCluster->globalX() - cluster->globalX());
+        if( abs(timeDifference) < 0.0001 ) correlationY[detectorID]->Fill(refCluster->globalY() - cluster->globalY());
         correlationTime[detectorID]->Fill( timeDifference );
         correlationTimeInt[detectorID]->Fill( timeDifferenceInt );
       }

branches/trunk/algorithms/Timepix3Clustering.C

@@ -22,15 +22,16 @@ int Timepix3Clustering::run(Clipboard* clipboard){
     // Check if they are a Timepix3
     string detectorID = parameters->detectors[det];
     if(parameters->detector[detectorID]->type() != "Timepix3") continue;
-		
+    tcout<<"Looking at detID "<<detectorID<<endl;
     // Get the pixels
     Timepix3Pixels* pixels = (Timepix3Pixels*)clipboard->get(detectorID,"pixels");
     if(pixels == NULL){
       tcout<<"Detector "<<detectorID<<" does not have any pixels on the clipboard"<<endl;
-      return 1;
+      continue;
     }
-    if(debug) tcout<<"Picked up "<<pixels->size()<<" pixels for device "<<detectorID<<endl;
-
+//    if(debug)
+      tcout<<"Picked up "<<pixels->size()<<" pixels for device "<<detectorID<<endl;
+    
     // Make the cluster storage
     Timepix3Clusters* deviceClusters = new Timepix3Clusters();
    
@@ -80,7 +81,7 @@ int Timepix3Clustering::run(Clipboard* clipboard){
     }
     
     // Put the clusters on the clipboard
-    clipboard->put(detectorID,"clusters",(TestBeamObjects*)deviceClusters);
+    if(deviceClusters->size() > 0) clipboard->put(detectorID,"clusters",(TestBeamObjects*)deviceClusters);
     if(debug) tcout<<"Made "<<deviceClusters->size()<<" clusters for device "<<detectorID<<endl;
 
   }
@@ -95,8 +96,8 @@ bool Timepix3Clustering::touching(Timepix3Pixel* neighbour,Timepix3Cluster* clus
   Timepix3Pixels pixels = cluster->pixels();
   for(int iPix=0;iPix<pixels.size();iPix++){
     
-    if( (pixels[iPix]->m_row - neighbour->m_row) <= 1 &&
-       (pixels[iPix]->m_column - neighbour->m_column) <= 1 ) Touching = true;
+    if( abs(pixels[iPix]->m_row - neighbour->m_row) <= 1 &&
+       abs(pixels[iPix]->m_column - neighbour->m_column) <= 1 ) Touching = true;
     
   }
   return Touching;

branches/trunk/algorithms/Timepix3EventLoader.C

@@ -83,6 +83,7 @@ void Timepix3EventLoader::initialise(Parameters* par){
 
 int Timepix3EventLoader::run(Clipboard* clipboard){
   
+  int endOfFiles = 0; int devices = 0;
 //  tcout<<"Current time is "<<parameters->currentTime<<endl;
   bool loadedData = false;
   // Loop through all registered detectors
@@ -101,11 +102,17 @@ int Timepix3EventLoader::run(Clipboard* clipboard){
 //      tcout<<"Loaded "<<deviceData->size()<<" pixels for device "<<detectorID<<endl;
       clipboard->put(detectorID,"pixels",(TestBeamObjects*)deviceData);
     }
+    
+    // Check if all devices have reached the end of file
+    devices++;
+    if(feof(m_currentFile[detectorID])) endOfFiles++;
   }
   
   // Increment the event time
   parameters->currentTime += parameters->eventLength;
   
+  if(endOfFiles == devices) return 0;
+  
   // If no data was loaded, tell the event loop to stop
   if(!loadedData) return 2;
   
@@ -141,8 +148,23 @@ 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 ){
 
+//  tcout<<"Loading data for device "<<detectorID<<endl;
+
   // Check if current file is open
-  if(m_currentFile[detectorID] == NULL){
+//  if(m_currentFile[detectorID] == NULL){
+  
+  if(m_currentFile[detectorID] == NULL || feof(m_currentFile[detectorID])){
+//    tcout<<"No current file open "<<endl;
+
+    
+    // If all files are finished, return
+    if(m_fileNumber[detectorID] == m_datafiles[detectorID].size()){
+//      tcout<<"All files have been analysed. There were "<<m_datafiles[detectorID].size()<<endl;
+      return false;
+    }
+    
+//    tcout<<"Opening file "<<m_fileNumber[detectorID]<<endl;
+    
     // Open a new file
     m_currentFile[detectorID] = fopen(m_datafiles[detectorID][m_fileNumber[detectorID]].c_str(),"rb");
     // Mark that this file is done
@@ -167,6 +189,10 @@ bool Timepix3EventLoader::loadData(string detectorID, Timepix3Pixels* devicedata
   // Now read the data packets.
   ULong64_t pixdata = 0; UShort_t thr = 0;
   int npixels=0;
+  bool fileNotFinished = false;
+  
+//  tcout<<"Ready to read data"<<endl;
+  
   // Read till the end of file (or till break)
   while (!feof(m_currentFile[detectorID])) {
     const int retval = fread(&pixdata, sizeof(ULong64_t), 1, m_currentFile[detectorID]);
@@ -210,6 +236,7 @@ bool Timepix3EventLoader::loadData(string detectorID, Timepix3Pixels* devicedata
       if( parameters->eventLength != 0. &&
          ((double)time/(4096. * 40000000.)) > (parameters->currentTime + parameters->eventLength) ){
         fseek(m_currentFile[detectorID], -1 * sizeof(ULong64_t), SEEK_CUR);
+        fileNotFinished = true;
         break;
       }
       
@@ -220,10 +247,18 @@ bool Timepix3EventLoader::loadData(string detectorID, Timepix3Pixels* devicedata
 
     }
     // Stop when we reach some large number of pixels (if events not based on time)
-    if(parameters->eventLength == 0. && npixels == 2000) break;
+    if(parameters->eventLength == 0. && npixels == 2000){
+      fileNotFinished = true;
+      break;
+    }
     
   }
   
+//  if(feof(m_currentFile[detectorID])){
+//    fclose(m_currentFile[detectorID]);
+//    m_currentFile[detectorID] == NULL;
+//    tcout<<"Closing file "<<endl;
+//  }
   if(npixels == 0) return false;
   
   return true;

branches/trunk/core/Analysis.C

@@ -27,9 +27,10 @@ void Analysis::run(){
   
   // Loop over all events, running each algorithm on each "event"
   cout << endl << "========================| Event loop |========================" << endl;
-  m_events=0;
+  m_events=1;
   while(1){
     bool run = true;
+    bool noData = false;
     cout<<"[Analysis] Running over event "<<m_events<<endl;
   	// Run all algorithms
     for(int algorithmNumber = 0; algorithmNumber<m_algorithms.size();algorithmNumber++) {
@@ -39,7 +40,7 @@ void Analysis::run(){
       m_algorithms[algorithmNumber]->getStopwatch()->Start(false);
       int check = m_algorithms[algorithmNumber]->run(m_clipboard);
       m_algorithms[algorithmNumber]->getStopwatch()->Stop();
-      if(check == 2) break; // Nothing to be done in this event
+      if(check == 2){noData = true; break;}// Nothing to be done in this event
       if(check == 0) run = false;
     }
     // Clear objects from this iteration from the clipboard
@@ -49,7 +50,7 @@ void Analysis::run(){
     // Check if we have reached the maximum number of events
     if(m_parameters->nEvents > 0 && m_events == m_parameters->nEvents) break;
     // Increment event number
-    m_events++;
+    if(!noData) m_events++;
   }
   
   // If running the gui, don't close until the user types a command

branches/trunk/core/Clipboard.h

@@ -44,6 +44,7 @@ public:
       int nElements = collection->size();
       for(TestBeamObjects::iterator it=collection->begin(); it!=collection->end(); it++) delete (*it);
       delete m_data[m_dataID[i]];
+      m_data.erase(m_dataID[i]);
     }
     m_dataID.clear();
   }