GRIFFINCollaboration / GRSISort

A lean, mean, sorting machine.
MIT License
21 stars 54 forks source link

No Bgo suppression applied to event-mixed coincidence for angular correlation matrices #1283

Closed Spagnole closed 8 months ago

Spagnole commented 2 years ago

Describe the bug I have been trying to perform angular correlation measurements with the 114Sn data (S1916). I have observed a staggering problem similar to what Fatima observed with the 80Ge data. While investigating this problem, I discover that for the event-mixed gamma-gamma coincidences, bgo suppression is not being applied. I first observed this by comparing the projections of gamma-gamma coincidence matrices between the event mixed and time-random coincidences. Shown below. EventMixing_NoBgo

To get a clearer picture of this issue I produced a simple selector file to compare singles spectra.

      void TestEventMixingSelector::FillHistograms() {
          //fill singles spectrum with addback and bgo suppression 
          for(auto g1 = 0; g1 < fGrif->GetSuppressedAddbackMultiplicity(fBgo); ++g1) {
              auto grif1 = fGrif->GetSuppressedAddbackHit(g1);
              if(PileUpRej) if(grif1->GetKValue() != kValue_grif) continue;
              fH2.at("gEA")->Fill(grif1->GetArrayNumber(), grif1->GetEnergy()); 
          }
          //fill separate spectrum using previous entry (mixed-event)
          for(unsigned int i = 0; i < grifQueue.size(); ++i) {
              fLastGrif = grifQueue[i];
              //fLastGrifBgo = grifBgoQueue[i];
              for(auto g2 = 0; g2 < fLastGrif.GetSuppressedAddbackMultiplicity(fBgo); ++g2) {
                  auto   grif2 = fLastGrif.GetSuppressedAddbackHit(g2);
                  if(PileUpRej) if(grif2->GetKValue() != kValue_grif) continue;
                  fH2.at("gEA_EM")->Fill(grif2->GetArrayNumber(), grif2->GetEnergy());
              }
          }     
          grifQueue.push_back(*fGrif); // add current entry to mixed-event list
          grifQueue.pop_front(); //remove previous entry to mixed-event list    
      }

This selector should fill a matrix of suppressed addback events from the TGriffin Branch (fGrif) of the current event in the tree. Then fill an identical matrix using the TGriffin Branch (fLastGrif) of the previous event in the tree. Then add fGrif to the list and remove the previous entry. The end result should be two identical matrices but this is not what I observe (see below). I have tried running the same script with GetAddbackHit (so no bgo suppression) and the resulting matrices are the same.

I think the suppression is not working as the fBgo events are not correlated with the fLastGrif events so the conditions for bgo suppression are never met.

EventMixing_NoBgo

I have tried to work around this issue by implementing a similar structure of past bgo events in the same way as the griffin events with some success.

void TestEventMixingSelector::FillHistograms() {

    for(auto g1 = 0; g1 < fGrif->GetSuppressedAddbackMultiplicity(fBgo); ++g1) {
        auto grif1 = fGrif->GetSuppressedAddbackHit(g1);
        if(PileUpRej) if(grif1->GetKValue() != kValue_grif) continue;
        fH2.at("gEA")->Fill(grif1->GetArrayNumber(), grif1->GetEnergy());
        grifHitQueue.push_back(*grif1);
    }
    for(unsigned int i = 0; i < grifQueue.size(); ++i) {
        fLastGrif = grifQueue[i];
        fLastGrifBgo = grifBgoQueue[i];
        for(auto g2 = 0; g2 < fLastGrif.GetSuppressedAddbackMultiplicity(fBgo); ++g2) {
            auto   grif2 = fLastGrif.GetSuppressedAddbackHit(g2);
            if(PileUpRej) if(grif2->GetKValue() != kValue_grif) continue;
            fH2.at("gEA_EM")->Fill(grif2->GetArrayNumber(), grif2->GetEnergy());
        }
        //work around for bgo suprression
        for(auto g2 = 0; g2 < fLastGrif.GetAddbackMultiplicity(); ++g2) {
            auto   grif2 = fLastGrif.GetSuppressedAddbackHit(g2);
            if(PileUpRej) if(grif2->GetKValue() != kValue_grif) continue;   
            bool bgo_flag = false;
            for(auto b1 = 0; b1 < fLastGrifBgo.GetMultiplicity(); b1++){
                auto bgo = fLastGrifBgo.GetBgoHit(b1);
                if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;
            }
            if( bgo_flag == true ) continue;                                
            fH2.at("gEA_EM2")->Fill(grif2->GetArrayNumber(), grif2->GetEnergy());
        }
    }   
    grifQueue.push_back(*fGrif);
    grifQueue.pop_front();
    grifBgoQueue.push_back(*fBgo);
    grifBgoQueue.pop_front();
}

This work around manages to reproduce the suppressed addback matrix but only if I limit the number of events in the tree. If I try to run the selector upon the whole tree, the process appears to crash on the matrices are empty. (see attached log)

TestEventMixingSelector19367_015.log

VinzenzBildstein commented 2 years ago

If you want to do suppressed event mixing you definitely need to keep the Griffin and BGO events together as you were doing in the last example. As for the crash, this line:

TGriffin Suppressed addback hits are out of range: vector::_M_range_check: __n (which is 0) >= this->size() (which is 0)
17:01:06  9893 Wrk-0.0 | *** Break ***: segmentation violation

says you are trying to access a suppressed add back hit that doesn't exist (in this case there are no suppressed addback hits. The issue is that one of your loops (under the comment "work around for bgo suprression") uses GetAddbackMultiplicity() and not GetSuppressedAddbackMultiplicity()!

VinzenzBildstein commented 2 years ago

Or alternatively that you use GetSuppressedAddbackHit(g2) instead of GetAddbackHit(g2) in that loop, depending on what you are trying to do there. And since you seem to do a manual suppression there you might want to use the unsuppressed hits for this?

Spagnole commented 2 years ago

Thanks for spotting the error for the GetAddbackMultiplicity and the GetSuppressedAddbackHit. I fixed that error and the script worked without errors.

To do the event-mixing for the angular correlation matrices I implemented the bgo suppression as I tried above

// Event mixing
for(unsigned int i = 0; i < grifQueue.size(); ++i) {
    fLastGrif = grifQueue[i];
    fLastGrifBgo = grifBgoQueue[i];
    for(auto g2 = 0; g2 < fLastGrif.GetAddbackMultiplicity(); ++g2) {
        //if(g1 == g2) continue;
        auto   grif2 = fLastGrif.GetAddbackHit(g2);
        if(PileUp) if(grif2->GetKValue() != KValueGrif) continue; //ignore pile up events
        if(CutCycle) if(fmod(grif2->GetTime()/1e9, MyCycleLength) < Cut) continue; //previous cutting condition
        if( !BeamOff(grif2) ) continue; //only use events in beam off window
        if( grif2->GetDetector() == 9 ) continue; //this is to remove problematic detector for testing purposes
        //flag for bgo suppression
        bool bgo_flag = false; 
        for(auto b1 = 0; b1 < fLastGrifBgo.GetMultiplicity(); b1++){
            auto bgo = fLastGrifBgo.GetBgoHit(b1);
            if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;
        }
        if( bgo_flag == true ) continue; //skipping events which fail bgo suppression   
        double angle = grif1->GetPosition(145.).Angle(grif2->GetPosition(145.)) * 180. / TMath::Pi();
        if(angle < 0.0001) continue;
        auto   angleIndex = fAngleMapAddback.lower_bound(angle - 0.0005);
        fH2["addbackAddbackHPMixed"]->Fill(grif1->GetArrayNumber(), grif2->GetArrayNumber());
        fH2["addbackAddbackMixed"]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        //reverse filling to make matrices symmetric
        fH2["addbackAddbackHPMixed"]->Fill(grif2->GetArrayNumber(), grif1->GetArrayNumber());
        fH2["addbackAddbackMixed"]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
        fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
        }
    }
}

grifQueue.push_back(*fGrif);
grifQueue.pop_front();
grifBgoQueue.push_back(*fBgo);
grifBgoQueue.pop_front();

This seems to work well. Looking at the comparison between the time-random coincidence matrix and the event-mixed matrix.

BgoSuppression_EventMixed

VinzenzBildstein commented 2 years ago

Does this now also solved the staggering issue?

Spagnole commented 2 years ago

I have not resorted all of the data yet after fixing this bgo suppression, it is a lot of data and will take some time.
However, I spoke with Fatima who observed staggering with the 80Ge data and she did not employ addback or bgo suppression.

Spagnole commented 2 years ago

The staggering issue looks to be resolved in this case after I removed one clover from the analysis (Grif 9 in this case). I removed this clover as its efficiency was around half compared to the other clovers.

Given that my issue was around the failure to perform bgo suppression for the event-mixed matrices and that issue seems to be resolved, I will close this issue.

Screenshot from 2022-08-15 16-38-42

magdasat commented 2 years ago

Hello!

I am now using Pietro's version of the angular correlation selector for my analysis and when I run it I get this error message:

[/net/npr2/d1/S1626/analysis_trees/analysis11502_000.root: adding ppg
Currently Running: AngularCorrelationSelectorAddback_2.C

Info in <TProofLite::SetQueryRunning>: starting query: 1
Info in <TProofQueryResult::SetRunning>: nwrks: 2
Info in <ACLiC>: unmodified script has already been compiled and loaded
Looking up for exact location of files: OK (1 files)         
Looking up for exact location of files: OK (1 files)         
Info in <TPacketizer::TPacketizer>: Initial number of workers: 2
Validating files: OK (1 files)         
0.0: caught exception triggered by signal '1' while processing dset:'TDSet:AnalysisTree', file:'/net/npr2/d1/S1626/analysis_trees/analysis11502_000.root' - check logs for possible stacktrace - last event: 2
Info in <TProofLite::MarkBad>: 
 +++ Message from master at nsl80.ph.liv.ac.uk : marking nsl80:-1 (0.0) as bad
 +++ Reason: undefined message in TProof::CollectInputFrom(...)

 +++ Message from master at nsl80.ph.liv.ac.uk : marking nsl80:-1 (0.0) as bad
 +++ Reason: undefined message in TProof::CollectInputFrom(...)

 +++ Most likely your code crashed
 +++ Please check the session logs for error messages either using
 +++ the 'Show logs' button or executing
 +++
 +++ root [] TProof::Mgr("nsl80.ph.liv.ac.uk")->GetSessionLogs()->Display("*")

0.1: caught exception triggered by signal '1' while processing dset:'TDSet:AnalysisTree', file:'/net/npr2/d1/S1626/analysis_trees/analysis11502_000.root' - check logs for possible stacktrace - last event: 82222
Info in <TProofLite::MarkBad>: 
 +++ Message from master at nsl80.ph.liv.ac.uk : marking nsl80:-1 (0.1) as bad
 +++ Reason: undefined message in TProof::CollectInputFrom(...)

 +++ Message from master at nsl80.ph.liv.ac.uk : marking nsl80:-1 (0.1) as bad
 +++ Reason: undefined message in TProof::CollectInputFrom(...)

 +++ Most likely your code crashed
 +++ Please check the session logs for error messages either using
 +++ the 'Show logs' button or executing
 +++
 +++ root [] TProof::Mgr("nsl80.ph.liv.ac.uk")->GetSessionLogs()->Display("*")

Using run 11502 subrun 0 3288060 events |>...................| 0.00 %
Writing TAnalysisOptions to New_AngularCorrelationAddback11502_000.root
Lite-0: all output objects have been merged                             

 *** Break *** segmentation violation
munmap_chunk(): invalid pointer
Aborted
mms@nsl80:/net/npr2/d1/S1626/mms/mms$]

At first I thought it was a memory issue and then that the number of the event-mixing matrices that I was defining in the .h file was too large. I made changes in both cases but the error message I was getting was the same. I then started commenting out parts of the code to see what is causing the issue and it appears that in the event-mixing section and inside the "for" loop of the bgo flag that Pietro has defined, this line

if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;

is confusing it for some reason. I am currently using version 3 of grsisort and I was wondering if the error can be related to the way that this mathematical expression is defined. It works fine for version 4 so maybe there is a different notation I should use? Also, the GetDetector() and GetTime() functions work fine outside this loop. This is the loop that I am referring to:

//flag for bgo suppression
            bool bgo_flag = false; 
            for(auto b1 = 0; b1 < fLastGrifBgo.GetMultiplicity(); b1++){
                auto bgo = fLastGrifBgo.GetBgoHit(b1);
                if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;
            }
            if( bgo_flag == true ) continue; //skipping events which fail bgo suppression   

Any ideas or comments you have are very much appriciated.

Thank you!!

VinzenzBildstein commented 2 years ago

I don't really know version 3 anymore, but if you narrowed it down to that part of the code, you can try and either put cout statements in to see what bgo is set too (maybe it's a null pointer?).

You can also check the log-files. There should be either a grsisort.log created in the directory you run grsiproof in, or it will tell you at the very end that there are log-files in ~/.proof//session--numbers/worker-0.0.log. You can also use the command find ~/.proof/ -name worker-0.0.log -exec ls -l {} \; to list all the log files you have and look at the most recent one.

Those log files might have more information on why the code crashes on that line.

As an aside: Why are you re-creating the BGO suppression yourself instead of using a loop over TGriffin::GetSuppressedMultiplicity() and using the hit TGriffin::GetSuppressedHit(<index>)?

VinzenzBildstein commented 2 years ago

Sorry @magdasat, I didn't see you question on slack. I'm not sure if there really is an issue with the suppression when using it correctly. The issue of the first code snippet that @Spagnole posted was that he tried suppressing fLastGrif with fBgo when it should have been fLastGrifBgo instead.

To understand this you need to realize that while we store GRIFFIN and BGO hits in two separate classes, they should never be split up. You always need to suppress the GRIFFIN detector with the BGO hits from the same event, you do not want to do event mixing here.

I'm not sure if @Spagnole tried going back to using the built-in suppression after fixing the other bugs in his selector?

Spagnole commented 2 years ago

So I have been building the angular correlation matrices to include bgo suppression in the event-mixing using the attached code. Firstly the event-mixed griffin and bgo events were declared in the header file as.

TGriffin fLastGrif;
std::deque<TGriffin> grifQueue;
TGriffinBgo fLastGrifBgo;
std::deque<TGriffinBgo> grifBgoQueue;

And the matrices were filled using the following

for(unsigned int i = 0; i < grifQueue.size(); ++i) {
    fLastGrif = grifQueue[i];
    fLastGrifBgo = grifBgoQueue[i];
    for(auto g2 = 0; g2 < fLastGrif.GetAddbackMultiplicity(); ++g2) {
        //if(g1 == g2) continue;
        auto   grif2 = fLastGrif.GetAddbackHit(g2);
        if(PileUp) if(grif2->GetKValue() != KValueGrif) continue; //ignore pile up events
        //if(CutCycle) if(fmod(grif2->GetTime()/1e9, MyCycleLength) < Cut) continue; //previous cutting condition
        if( !InCyle(grif2) ) continue; //only use events in beam off window
        //if( grif2->GetDetector() == 9 ) continue; //this is to remove problematic detector for testing purposes
        //flag for bgo suppression
        bool bgo_flag = false; 
        for(auto b1 = 0; b1 < fLastGrifBgo.GetMultiplicity(); b1++){
            auto bgo = fLastGrifBgo.GetBgoHit(b1);
            if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;
        }
        if( bgo_flag == true ) continue; //skipping events which fail bgo suppression   
        double angle = grif1->GetPosition(145.).Angle(grif2->GetPosition(145.)) * 180. / TMath::Pi();
        if(angle < 0.0001) continue;
        auto   angleIndex = fAngleMapAddback.lower_bound(angle - 0.0005);
        fH2["addbackAddbackHPMixed"]->Fill(grif1->GetArrayNumber(), grif2->GetArrayNumber());
        fH2["addbackAddbackMixed"]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        //reverse filling to make matrices symmetric
        //fH2["addbackAddbackHPMixed"]->Fill(grif2->GetArrayNumber(), grif1->GetArrayNumber());
        //fH2["addbackAddbackMixed"]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
        //fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
    }
}

grifQueue.push_back(*fGrif);
grifQueue.pop_front();
grifBgoQueue.push_back(*fBgo);
grifBgoQueue.pop_front();

In order to perform bgo suppression with the same method as for the current events (fGrif and fBgo), I had to change fLastGrif and grifBgoQueue to pointers.

TGriffin *fLastGrif;
std::deque<TGriffin*> grifQueue;
TGriffinBgo *fLastGrifBgo;
std::deque<TGriffinBgo*> grifBgoQueue;

Then the bgo suppression is done as follows.

for(unsigned int i = 0; i < grifQueue.size(); ++i) {
    fLastGrif = grifQueue[i];
    fLastGrifBgo = grifBgoQueue[i];
    for(auto g2 = 0; g2 < fLastGrif->GetSuppressedAddbackMultiplicity(fLastGrifBgo); ++g2) {
        //if(g1 == g2) continue;
        auto   grif2 = fLastGrif->GetSuppressedAddbackHit(g2);
        if(PileUp) if(grif2->GetKValue() != KValueGrif) continue; //ignore pile up events
        //if(CutCycle) if(fmod(grif2->GetTime()/1e9, MyCycleLength) < Cut) continue; //previous cutting condition
        if( !InCyle(grif2) ) continue; //only use events in beam off window
        //if( grif2->GetDetector() == 9 ) continue; //this is to remove problematic detector for testing purposes
        double angle = grif1->GetPosition(145.).Angle(grif2->GetPosition(145.)) * 180. / TMath::Pi();
        if(angle < 0.0001) continue;
        auto   angleIndex = fAngleMapAddback.lower_bound(angle - 0.0005);
        fH2["addbackAddbackHPMixed"]->Fill(grif1->GetArrayNumber(), grif2->GetArrayNumber());
        fH2["addbackAddbackMixed"]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
        //reverse filling to make matrices symmetric
        //fH2["addbackAddbackHPMixed"]->Fill(grif2->GetArrayNumber(), grif1->GetArrayNumber());
        //fH2["addbackAddbackMixed"]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
        //fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
        }
    }
}

grifQueue.push_back(fGrif);
grifQueue.pop_front();
grifBgoQueue.push_back(fBgo);
grifBgoQueue.pop_front();

I have attached a comparison of the results below. I have sorted the same file with both version of the code and sorted the file again with the first version of the code but I removed the bgo suppression by commenting out the line

if( bgo_flag == true ) continue; //skipping events which fail bgo suppression

Additionally the matrices are not being filled symmetrically for this comparison. This is so that the I can just day the projection of the event-mixed energies.

The histogram produced from the first code is drawn in blue, the second code in red and with no bgo suppression applied is drawn in magenta. I am not sure why there is such is difference in the results between the first and second code.

Comparrison

I have been building the angular correlation matrices with the first code. In the previous angular correlation matrices I produced with the second code, it appeared that no bgo suppression was being applied. I am not sure what has changed to fix this. I have just written the second code fresh so perhaps I made a mistake in my previous selector.

mrocchini commented 2 years ago

Hi Pietro. Would it be possible to fit a few lines in the spectra "first code" and "second code" to see if we are losing good data or if we simply have a better suppression with the second code? I see you have some pretty strong lines at ~500 keV, ~800 keV, 1400 keV, 1700 keV and 2800 keV. In the former case (losing good data), it would be interesting to see also how much we are losing at the different energies. Thanks.

Spagnole commented 2 years ago

Here is a table of peak intensities comparing the first and second code with no bgo suppression applied. For the peak intensities with the first code, they are very similar to the intensities without bgo suppression with the exception of the 511-keV peak. This is not the case for the second code.

<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">

  | No Bgo Suppression |   | First Code |   |   |   | Second Code |   |   |   -- | -- | -- | -- | -- | -- | -- | -- | -- | -- | -- Energy | Int | err | Int | err | Int (First code / No Supp) | err | Int | err | Int (Second code / No Supp) | err 393.403 | 53690 | 1412 | 51582 | 463 | 0.961 | 0.027 | 54940 | 2595 | 1.023 | 0.055 510.981 | 173848 | 744 | 108157 | 567 | 0.622 | 0.004 | 90691 | 913 | 0.522 | 0.006 569.835 | 217553 | 609 | 208839 | 542 | 0.960 | 0.004 | 227176 | 550 | 1.044 | 0.004 814.962 | 1098070 | 1080 | 1069910 | 1113 | 0.974 | 0.001 | 799868 | 914 | 0.728 | 0.001 963.426 | 50235 | 347 | 48848 | 324 | 0.972 | 0.009 | 45395 | 2113 | 0.904 | 0.043 1384.27 | 284265 | 605 | 282047 | 2987 | 0.992 | 0.011 | 130908 | 407 | 0.461 | 0.002 1712.29 | 93121 | 379 | 91256 | 1724 | 0.980 | 0.019 | 68643 | 294 | 0.737 | 0.004 2821.03 | 101083 | 352 | 99934 | 330 | 0.989 | 0.005 | 45016 | 223 | 0.445 | 0.003
Spagnole commented 2 years ago

@magdasat I believe I managed to get the bgo suppression to work with grsisort version 3. I have attached a figure below comparing with and without bgo suppression. Again to ignore bgo suppression I commented out the flag. if( bgo_flag == true ) continue; //skipping events which fail bgo suppression

I used a similar method as my "First Code" discussed above. A copy of the selector file is provided below.

Firstly I had to convert fLastGrif and fLastGrifBgo to pointers

TGriffin *fLastGrif;
  std::deque<TGriffin*> grifQueue;
  TGriffinBgo *fLastGrifBgo;
  std::deque<TGriffinBgo*> grifBgoQueue;
for(unsigned int i = 0; i < grifQueue.size(); ++i) {
fLastGrif = grifQueue[i];
fLastGrifBgo = grifBgoQueue[i];
for(auto g2 = 0; g2 < fLastGrif->GetAddbackMultiplicity(); ++g2) {
    //if(g1 == g2) continue;
    auto   grif2 = fLastGrif->GetAddbackHit(g2);
    //flag for bgo suppression
    bool bgo_flag = false; 
    for(auto b1 = 0; b1 < fLastGrifBgo->GetMultiplicity(); b1++){
        auto bgo = fLastGrifBgo->GetBgoHit(b1);
        //cout << bgo->GetDetector() << endl;
        fH1["BgoEvents"]->Fill(bgo->GetDetector()); 
        fH2["BgoGrifMap"]->Fill(grif2->GetDetector(), bgo->GetDetector()); 
        if( bgo->GetDetector() == grif2->GetDetector() ){ 
            double td =  bgo->GetTime() - grif2->GetTime();
            fH1["BgoGrif_td"]->Fill(td);    
        }   
        if( bgo->GetDetector() == grif2->GetDetector() && TMath::Abs(bgo->GetTime() - grif2->GetTime() ) < 300. ) bgo_flag = true;
    }
    //if( bgo_flag == true ) continue; //skipping events which fail bgo suppression 
    double angle = grif1->GetPosition(145.).Angle(grif2->GetPosition(145.)) * 180. / TMath::Pi();
    if(angle < 0.0001) continue;
    auto   angleIndex = fAngleMapAddback.lower_bound(angle - 0.0005);
    fH2["addbackAddbackHPMixed"]->Fill(grif1->GetArrayNumber(), grif2->GetArrayNumber());
    fH2["addbackAddbackMixed"]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
    fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif1->GetEnergy(), grif2->GetEnergy());
    //reverse filling to make matrices symmetric
    //fH2["addbackAddbackHPMixed"]->Fill(grif2->GetArrayNumber(), grif1->GetArrayNumber());
    //fH2["addbackAddbackMixed"]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
    //fH2[Form("addbackAddbackMixed%d", angleIndex->second)]->Fill(grif2->GetEnergy(), grif1->GetEnergy());
}

}

grifQueue.push_back(fGrif);
grifQueue.pop_front();
grifBgoQueue.push_back(fBgo);
grifBgoQueue.pop_front();

EventMixingBgoSuppV3

MyNewAngularCorrelationSelector.h.txt MyNewAngularCorrelationSelector.C.txt

VinzenzBildstein commented 2 years ago

Just a general comment here, using the std::deque with objects (instead of pointers) definitely introduces a memory leak. This is because removing an element from the deque does not guarantee that the destructor of that element is called.

mrocchini commented 2 years ago

Hi Pietro, this is interesting. Thanks for checking. I will try to cross-check with other data. I should be able to work on this in the next month or so.

VinzenzBildstein commented 1 year ago

I'm assuming this issue has been solved and we can close it? If needed we can re-open it.

magdasat commented 1 year ago

Hi Vinzenz! yes, you can close it thank you. I am still trying things so I will comment here again as soon as I have an update.

lpgaff commented 9 months ago

Just a general comment here, using the std::deque with objects (instead of pointers) definitely introduces a memory leak. This is because removing an element from the deque does not guarantee that the destructor of that element is called.

Sorry to drag up an old post, but we have an issue with the std::deque for the TGriffin and TGriffinBgo events.

SYSTEM = Linux
OS = Debian GNU/Linux
VER = 11
ROOT VERSION = 6.28/04
COMMIT = [09e52e9](https://github.com/GRIFFINCollaboration/GRSISort/commit/09e52e93a4dc453bec3bbf042156afd243446ef2)

Vinzenz is of course correct that using the std::deque with pointers gives a memory leak because they are not destroyed with pop_front. However, I am still experiencing memory leaks without using pointers. I've tried so many crazy ways to get around this, but the memory just keeps increasing until grsiproof crashes if we use the push_back and pop_front combination of commands, including doing a "manual" version of pushing and poping. If I comment this part out, the memory usage remains constant throughout the execution of the code.

This is how it is done currently:

grifQueue.push_back(*fGrif);
grifQueue.pop_front();
grifBgoQueue.push_back(*fBgo);
grifBgoQueue.pop_front();

I suspect it has its origins in one of two places:

The latter would probably rear it's head elsewhere, so I am leaning towards one of the first two. Any suggestions? Or can anybody reproduce this beahviour by watching the memory usage in top or similar during execution of grsiproof on this selector?

P.S. Reducing the Nmixed parameter sufficiently can allow us to finish a file before hitting the memory limit of the PC, but that isn't really a satisfactory solution to a memory leak :-)

VinzenzBildstein commented 9 months ago

No problem dragging up an old issue, but next time please also re-open it.

The memory leak you are seeing is exactly the kind of issue I mentioned. Looking at the first two lines, what you are doing is:

You are assuming that removing the copy from the deque also frees the memory, but it does not! That's where your memory leak is coming from.

So if you want to use a deque, do something like this:

   grifQueue.emplace_back(new TGriffin(*fGrif));
   delete grifQueue.front();
   grifQueue.pop_front();

where grifQueue is std::deque<TGriffin*> fGrifQueue; so a deque of pointers, and not objects. I think an example of how you can do this is in the AngularCorrelationSelector.C file (though I notice that this one is missing the deque for the BGOs).

VinzenzBildstein commented 9 months ago

Just to clarify: you DO want to use a deque with pointers (where you can take care of deleting the memory properly). And you DON'T want to use a deque of objects (where you can not rely on their destructors being called).

lpgaff commented 9 months ago

Ahhaaaaaa. You are a hero, Vinzenz. Thank you. The missing thing was to delete the pointer at the front of the deque before popping.

Just FYI, the AngularCorrelationSelector in the main branch does not use a std::deque for event mixing, it just uses the previous event, i.e Nmixed = 1. See here: https://github.com/GRIFFINCollaboration/GRSISort/blob/main/GRSIProof/AngularCorrelationSelector.C

Maybe that should be updated (yes, I realise I just gave myself a task to do...)

VinzenzBildstein commented 9 months ago

Yes, the normal angular correlation selector (without add back) doesn't use multiple mixed events. Nor does the angular correlation helper. I'm also not sure if any of the other issues that have come up with the angular correlations have been fixed in the ones that are part of the repo.

VinzenzBildstein commented 8 months ago

I've updated the angular correlation helper. It now uses deques to store multiple TGriffin and TGriffinBgo from previous events. I will close this issue, feel free to re-open it if you find other issues with angular correlations.