Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update to StFcsWaveformFitMaker to use new algorithm #401

Merged
merged 11 commits into from
Nov 17, 2022

Conversation

dkapukchyan
Copy link
Contributor

This is an update to the StFcsWaveformFitMaker algorithm to improve the peak finding method and energy determination from the raw signal data.

}
if(mMeasureTime){
mTime=new TH1F("FitTime","FitTime; msec",200,0,6);
mTime=new TH1F("FitTime","FitTime;msec",1000,0,1000);
Copy link
Member

@plexoos plexoos Sep 19, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am using this line to comment on a common problem which this PR seems to introduce in quite a few other places.
The mTime histogram can be allocated multiple times during a job when the framework calls InitRun(), however, the corresponding object is deleted only once in the destructor, hence causing a memory leak. Depending on the logic you are trying to implement here, you may need to allocate histograms only once in the Init() method. It may be also easier to rely on the automatic object lifetime than trying to manage the memory manually like you do.

Copy link
Contributor Author

@dkapukchyan dkapukchyan Sep 20, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought InitRun was only called once during a job, which is why I put the histogram creations there. I agree that I only need these histograms created once per job. Where should I move the histogram allocations so that they only get created once? Hopefully I can move these before this pull request is made so that I don't need to make a new one.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Move to Init(). In principle, InitRun can be called multiple times per job if you are reading file(s) from multiple runs.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If it needs to be done only once, is there a reason not to initialize the histograms in the constructor? Also, if there is no good reason to allocate them with the new operator, you can avoid the many delete statements in the destructor.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

NB: During StMaker::Init(), any newly instantiated histograms are added to the maker's "hist" subdirectory (in TDataset terms) which allows them to be found using the TDataset methods and enables inclusion in the BFC's hist.root output files.

// Add the Maker histograms in the Maker histograms list
// and remove it from the ROOT system directory
gROOT->cd();
TIter nextHist(gDirectory->GetList());
Int_t ready = !objLast;
while((objHist=nextHist())) {// loop over gDirectory
if (!ready && objHist!=objLast) continue;
ready = 1999;
if (objHist==objLast) continue;
if (!objHist->InheritsFrom("TH1")) continue;

// Move the histogram from the ROOT list into the "maker's" list
((TH1*)objHist)->SetDirectory(0);
maker->AddHist((TH1*)objHist);
}

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Moving them to the constructor means that these histograms will always be created, whereas now they are only created when a certain test level is set using setTest(). I thought this would be best since in the production these histograms don't need to be created as they are mostly used to quantify how well the algorithm is working. If I move them to the constructor then I would have to add an argument for the test level and file name where those histograms will be saved. I think doing that will only make this class more complicated so I left them in Init(). That way whoever calls this class can turn testing on or off before it hits the chain.

Copy link
Member

@plexoos plexoos left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

StFcsWaveformFitMaker is already included in the fcs chain option. According to this page https://www.star.bnl.gov/devcgi/dbProdOptionRetrv.pl it was used once in
pp2022 StiCA BEmcChkStat epdhit -hitfilt for "Preliminary chain for run 2022 fastoffline production". From this I conclude that changing the algorithm should not have significant consequences

StRoot/StFcsWaveformFitMaker/PeakWindow.cxx Outdated Show resolved Hide resolved
StRoot/StFcsWaveformFitMaker/PeakWindow.cxx Outdated Show resolved Hide resolved
StRoot/StFcsWaveformFitMaker/PeakAnaVirtualPainter.h Outdated Show resolved Hide resolved
StRoot/StFcsWaveformFitMaker/PeakAnaPainter.h Outdated Show resolved Hide resolved
StRoot/StFcsDbMaker/StFcsDbMaker.h Show resolved Hide resolved
{
graph->SetPoint(ipeak, Ana.GetPeak(ipeak).mPeakX, Ana.GetPeak(ipeak).MidPoint() );
}
PeakAna* NewAna = new PeakAna(Ana,graph);
Copy link
Member

@plexoos plexoos Sep 27, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This can cause a memory leak. Why not return the object? EDIT: I mean "return by value"

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not sure whether you are referring to graph or NewAna. I don't think graph will give a memory leak since the next line calls NewAna->ForceInternal(). This means that the newly created graph object is "owned" by NewAna and it will delete that newly created graph when it is destroyed or SetData() is called. There is no way to change mInternalSignal to false after that so NewAna will never think that it doesn't own graph. If you are referring to NewAna then I agree it will if someone forgets to delete it. I don't want to return the graph object itself because if someone does PeakAna* newana = new PeakAna(OldAna,ConvertPeaksToAna(OldAna)) then the way this constructor works it will think the graph object is external to it and it won't delete it and that will be a memory leak. Which means that the caller will have to remember to always call newana->ForceInternal() in this situation. Here I called it automatically to prevent such a situation so that's why I don't return the graph object.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually, I meant you could return the NewAna object by value rather than by pointer, not the graph object. I am a bit lost on whether PeakAna owns the graph object or not and under what conditions... Hopefully, this complex (for me) logic is correct and there are no memory leaks.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought about returning by value but my understanding is that when you return by value it returns a copy of the object not the object created in the function itself. This is not something I want for a large object like PeakAna. Also, if you try to get around this by returning a reference instead then the returned reference is to the object created in the function which has now been invalidated and hence does not exist. This is why I return the pointer to the new object. Is my understanding incorrect and is there a safe way to do this without copying the object?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is safe to assume that your compiler implements the return value optimization. Moreover, copy elision is guaranteed when an object is returned directly (as of C++17)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is it going to be addressed?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

done

StRoot/StFcsWaveformFitMaker/StFcsWaveformFitMaker.h Outdated Show resolved Hide resolved
Comment on lines 75 to 79
virtual int Init();
virtual int InitRun(int runNumber);
virtual int Make();
virtual int Finish();
virtual void Clear(Option_t* option = "");
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What you should really use here is the override keyword rather than virtual

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I tried using override and compiling with cons and it complains that symbol override is not defined. Can you provide a working example?

StRoot/StFcsWaveformFitMaker/StFcsWaveformFitMaker.h Outdated Show resolved Hide resolved
StRoot/StFcsWaveformFitMaker/StFcsWaveformFitMaker.h Outdated Show resolved Hide resolved
@starsdong
Copy link
Member

Akio, Daniel, we should try to conclude on this PR. As codeowners, could you please review and approve this PR?
Dmitri, in the meantime, the PR seems to be out-of-sync with main branch. Could you update the branch?

@dkapukchyan
Copy link
Contributor Author

Xilin looked at invariant mass plot and finds no difference in the pi0 peak position and shape using the StFcsWaveformFitMaker in this commit and current StFcsWaveformFitMaker.

Copy link
Contributor

@akioogawa akioogawa left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is now well tested and ready to go.

Comment on lines 120 to 123
if(mH1_Baseline!=0) { delete mH1_Baseline; }
if(mF1_SignalFit!=0) { delete mF1_SignalFit; }
if(mF1_BaselineFit!=0) { delete mF1_BaselineFit; }

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
if(mH1_Baseline!=0) { delete mH1_Baseline; }
if(mF1_SignalFit!=0) { delete mF1_SignalFit; }
if(mF1_BaselineFit!=0) { delete mF1_BaselineFit; }
delete mH1_Baseline;
delete mF1_SignalFit;
delete mF1_BaselineFit;

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

done

Comment on lines 444 to 446
//if( peak>0 && radc<=last) peak ^= 1<<rightshift ;
//std::cout << " - SumDep0Mod::Case2|peak:"<<peak << "|radc:"<<radc<<"|last:"<<last << std::endl;
//if( peak!=0 || radc<last) peak = -1 ;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you please clean up the code and remove the commented code in all modified files? These distract focus from the important lines

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I removed what I could. Please let me know if there is more that I missed. Long comments describing what the algorithms do and how they work I left in since I don't know where else to put them.

{
graph->SetPoint(ipeak, Ana.GetPeak(ipeak).mPeakX, Ana.GetPeak(ipeak).MidPoint() );
}
PeakAna* NewAna = new PeakAna(Ana,graph);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is it going to be addressed?

Comment on lines +267 to +269
std::cout << "|P::Graph:"<<mG_Data;
std::cout << "|RangeX:["<<mXRangeMin<<","<<mXRangeMax<<"]";
std::cout << "|RangeY:["<<mYRangeMin<<","<<mYRangeMax<<"]";
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is it going to be addressed?


PeakAna* PeakAna::MeanFilter( Int_t sizeavgs, bool copy )
{
if( this->GetData()==0 ){return this;}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You use this-> a lot for class members. In most of our code it is usually omitted. I suggest you do the same.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Removed what I found. I only left the ones relating to virtual methods.

StRoot/StFcsWaveformFitMaker/StFcsPulseAna.h Show resolved Hide resolved
Comment on lines +418 to +424
if( mH1_NPeaksAkio!=0 ){ mH1_NPeaksAkio->Write(); }
if( mH1_NPeaksFilteredAkio!=0 ){ mH1_NPeaksFilteredAkio->Write(); }
if( mH1_PeakTiming!=0 ){mH1_PeakTiming->Write();}

if( mH2_NPeakvsPeakXdiff!=0 ) {mH2_NPeakvsPeakXdiff->Write();}
if( mH2_NPeakvsPeakYratio!=0) {mH2_NPeakvsPeakYratio->Write();}
if( mH1_VOverlap!=0 ) {mH1_VOverlap->Write();}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Have you tried to make a single call to rootFile->Write(); instead of doing this for individual objects?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have tried it before for other classes and it doesn't always work. It probably has to do with the order in which you create TFile and histograms so I explicitly write them in now. If it is really a big problem I can try to see if rootFile->Write() will work in this instance.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, please try again. ROOT has this notion of a current directory so, when you open a new file it switches to that current directory and it becomes a parent for all TObject-s created after that. The common trick with ROOT is to save the current directory and restore it later when you are done creating your histograms. Since you do it in a Maker this can interfere with other opened files, although I recall you have a switch which disables histogram filling by default. Since no tests were provided for all of this functionality we can only hope that it all works as intended.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I tried with just rootFile->Write() and it didn't work. I understand the issue but I don't know what else to do to avoid this. You are correct that if no test level is called then no output needs to be saved. The test levels are for single file tests and should never be called in production. They are intended to help understand and improve the algorithm rather than for generating generic QA histograms for StFcsWaveformFitMaker.

Comment on lines +418 to +424
if( mH1_NPeaksAkio!=0 ){ mH1_NPeaksAkio->Write(); }
if( mH1_NPeaksFilteredAkio!=0 ){ mH1_NPeaksFilteredAkio->Write(); }
if( mH1_PeakTiming!=0 ){mH1_PeakTiming->Write();}

if( mH2_NPeakvsPeakXdiff!=0 ) {mH2_NPeakvsPeakXdiff->Write();}
if( mH2_NPeakvsPeakYratio!=0) {mH2_NPeakvsPeakYratio->Write();}
if( mH1_VOverlap!=0 ) {mH1_VOverlap->Write();}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, please try again. ROOT has this notion of a current directory so, when you open a new file it switches to that current directory and it becomes a parent for all TObject-s created after that. The common trick with ROOT is to save the current directory and restore it later when you are done creating your histograms. Since you do it in a Maker this can interfere with other opened files, although I recall you have a switch which disables histogram filling by default. Since no tests were provided for all of this functionality we can only hope that it all works as intended.

Comment on lines +267 to +269
std::cout << "|P::Graph:"<<mG_Data;
std::cout << "|RangeX:["<<mXRangeMin<<","<<mXRangeMax<<"]";
std::cout << "|RangeY:["<<mYRangeMin<<","<<mYRangeMax<<"]";
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not aware of any exclusion to the rule requiring the use of LOG_ stream objects. The advantage is that we can turn off or limit logging to a predefined level. I hope in production only errors and maybe warnings are logged to avoid unnecessary writes to, for example, LOG_DEBUG and LOG_INFO streams

Comment on lines 181 to 183
if( mH1_PeakTimingGaus!=0 ) {delete mH1_PeakTimingGaus;}
if( mH1_PeakTimingPuls!=0 ) {delete mH1_PeakTimingPuls;}
if( mH2_PeakTimingCompare!=0 ) {delete mH2_PeakTimingCompare;}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought we already discussed this: there is no need to check for the pointer to be non null before deleting the object. Isn't nicer to have less code than more code?

Suggested change
if( mH1_PeakTimingGaus!=0 ) {delete mH1_PeakTimingGaus;}
if( mH1_PeakTimingPuls!=0 ) {delete mH1_PeakTimingPuls;}
if( mH2_PeakTimingCompare!=0 ) {delete mH2_PeakTimingCompare;}
delete mH1_PeakTimingGaus;
delete mH1_PeakTimingPuls;
delete mH2_PeakTimingCompare;

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Done. Sorry I had forgotten about this file :)

@jdbrice jdbrice merged commit 347894e into star-bnl:main Nov 17, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

7 participants