Thursday, May 30, 2019
NN problem with detector features
Writing features files from flat trees output by DSelector. Should be easy enough to keep detector information, BUT not every event has meaningful output for each of these. E.g., for events where the K+ does not enter the BCAL, BCAL values are either 0 or 1000 (the latter a flag value, I assume). A NN will likely have problems handling this (feature standardization will be screwy), but a BDT should be able to do it easily.
Tuesday, May 28, 2019
Next up: features files from flat trees
Next step is to write code that takes the flat trees output by the protonTRUTH DSelector and produces features files. Should give features that I didn't include in previous features files some type of tag that allows them to be removed from the dataframes prior to fitting.
Make sure to OBEY THE SAME FEATURE NAMES from previous scheme. Otherwise will have to change lots of NN code downstream.
A thing to worry about: these flat trees have no incorrect permutation events. Because proton/K+ separation is not good above a certain track momentum, this is unlikely to be true of actual data for which we can't cut on particle TRUTH id. Need to keep this in mind. Perhaps another NN that handles positive-track p/K+ id?
Make sure to OBEY THE SAME FEATURE NAMES from previous scheme. Otherwise will have to change lots of NN code downstream.
A thing to worry about: these flat trees have no incorrect permutation events. Because proton/K+ separation is not good above a certain track momentum, this is unlikely to be true of actual data for which we can't cut on particle TRUTH id. Need to keep this in mind. Perhaps another NN that handles positive-track p/K+ id?
NN training idea...
Possible to add result figure of merit (i.e. opposite of Punzi, significance) to the loss function when training neural networks? Since we know the approximate ratio of signal to background reactions, it is possible to predict S and B for a specific cut.
After more thought, this will be hard to do, as S and B will not update with every event processed (or batch of events), and thus will not contribute to back propagation.
After more thought, this will be hard to do, as S and B will not update with every event processed (or batch of events), and thus will not contribute to back propagation.
Updated fastpion hdgeant4
Rebuilt hdgeant4 pion zapper last night. Submitted jobs around 9:00, and they sat in queue for a while. Just finishing now.
Looking at the resulting vertex files, definitely more pions decaying in the detector, but still probably not enough. Scale the decay length down to 50 and retry...
UPDATE: reran with the max length for each detector element set to 50. This did not produce many more pi- decays. Might have to rethink how I'm handling this...
This is even stranger now that I'm remembering that the lengths in hdgeant4 are in mm...
Resetting max length to 100 mm, but changing max track id from 3 to 6:
if (track.GetTrackID() < 6 &&
Rerunning...
UPDATE: Yep, that did it! Pretty much all of the pions are picked up by this change. Resetting the decay length to 500 mm and rerunning...
UPDATE: Looking good!
NOW that that's all tidied up, will need a new in-depth look at the fastpi MC. This will likely have consequences for the ways that vertexing and NN classification were working out, since the fastpi MC wasn't really and different from the regular ppim MC. Will also have to have Megan histogram the pion decay locations.
Looking at the resulting vertex files, definitely more pions decaying in the detector, but still probably not enough. Scale the decay length down to 50 and retry...
UPDATE: reran with the max length for each detector element set to 50. This did not produce many more pi- decays. Might have to rethink how I'm handling this...
This is even stranger now that I'm remembering that the lengths in hdgeant4 are in mm...
Resetting max length to 100 mm, but changing max track id from 3 to 6:
if (track.GetTrackID() < 6 &&
Rerunning...
UPDATE: Yep, that did it! Pretty much all of the pions are picked up by this change. Resetting the decay length to 500 mm and rerunning...
UPDATE: Looking good!
NOW that that's all tidied up, will need a new in-depth look at the fastpi MC. This will likely have consequences for the ways that vertexing and NN classification were working out, since the fastpi MC wasn't really and different from the regular ppim MC. Will also have to have Megan histogram the pion decay locations.
Monday, May 27, 2019
fast pions not fast enough... Vertex files for Megan
Good news is that I think that the hddm to vertex file generator is working. Assembling files to send to Megan. Decided against putting the decaying particle names in the vertex files b/c it was too cumbersome to switch for the different topos AND it's inefficient. The particle id numbers in the file refer to the particle ids in the raw files.
In checking these files, it looks like only about 40% of the fast pion events have pions that decay (at all). This could be due to a) not having a short-enough lifetime or b) track.GetTrackID() < 3 not being met for the pions in MyPrimaryPionZapper.cc in hdgeant4.
Looking at the code, it is probably the former, as I currently have the path length set to be up to 5 m! Yikes:
return (rand() % 1000)/1000.0 * 500;
Turning this down to 1 m:
return (rand() % 1000)/1000.0 * 100;
... and rebuilding hdgeant4 and rerunning fastpi mc.
In checking these files, it looks like only about 40% of the fast pion events have pions that decay (at all). This could be due to a) not having a short-enough lifetime or b) track.GetTrackID() < 3 not being met for the pions in MyPrimaryPionZapper.cc in hdgeant4.
Looking at the code, it is probably the former, as I currently have the path length set to be up to 5 m! Yikes:
return (rand() % 1000)/1000.0 * 500;
Turning this down to 1 m:
return (rand() % 1000)/1000.0 * 100;
... and rebuilding hdgeant4 and rerunning fastpi mc.
Sunday, May 26, 2019
New problem with MMoffK+ spectrum difference between sig and bkgd [SOLVED]
See below. Significant difference between Lambda mass from MMoffK+ when comparing signal and background. Need to check raw MC generator.
Checking:
1. kL_phasesp_genr8.py. Correctly assigning particle type -18 to all Lambdas (so that GEANT doesn't track/decay the Lambda).
LIGHTNING BOLT -- this is probably OK, the result of the KF trying to make the background events look like sl_mu events. Check the pre-KF MMoffK+:
Yessss
Checking:
1. kL_phasesp_genr8.py. Correctly assigning particle type -18 to all Lambdas (so that GEANT doesn't track/decay the Lambda).
LIGHTNING BOLT -- this is probably OK, the result of the KF trying to make the background events look like sl_mu events. Check the pre-KF MMoffK+:
Yessss
DSelector TRUTH combos
For now, it looks like I should cut on the truth ID for the proton only. This will keep all three types of decays, AND it will keep events for which the K+ decays in the detector and is tracked like a mu+. As of right now, there doesn't seem to be a big difference between keeping and cutting these latter events, but will investigate further soon.
Attached is a plot of the MMoffK+ distribution pre- and post the selection of the proton thrown index.
Attached is a plot of the MMoffK+ distribution pre- and post the selection of the proton thrown index.
Saturday, May 25, 2019
DSelector flat trees
Flat trees from DSelector are working nicely. All of the quantities that we are likely to want are easy to access. Switching to flat trees from now on. Will write a piece of code that rips a flat tree to ascii features files for later analysis.
Friday, May 24, 2019
DSelector TRUTH combo output
Spent some time working on getting DSelector to output only combos that match the truth PIDs. As Sean suggested, IsTrueCombo() didn't work well... or at all. It always returned false. This might have something to do with how I set the particle indexes in my raw file. I set K+, p, and mu- to have idx 1, 2, 3; but these idx are 0, 1, 2 when accessed in DSelector. Might be a mis-match there.
I coded a workaround, which is to check each particle wrapper's ThrownID and compare to 1, 2, 3. In order to handle all three decay topologies, should be enough to check the positive particle permutations (K+, p). This way, the code will be agnostic about whether the tracked negative particle is a pi- or mu-.
After implementation, this seems to do the trick -- the combinatorical tail in the Lambda invariant mass is dramatically decreased.
I want to add the track thrown IDs to the output tree, but this is causing a crash right now. It might be better to go to a flat tree for the time being, as adding flat custom branches is less weird.
I coded a workaround, which is to check each particle wrapper's ThrownID and compare to 1, 2, 3. In order to handle all three decay topologies, should be enough to check the positive particle permutations (K+, p). This way, the code will be agnostic about whether the tracked negative particle is a pi- or mu-.
After implementation, this seems to do the trick -- the combinatorical tail in the Lambda invariant mass is dramatically decreased.
I want to add the track thrown IDs to the output tree, but this is causing a crash right now. It might be better to go to a flat tree for the time being, as adding flat custom branches is less weird.
Thursday, May 23, 2019
To understand next...
500 files of each MC flavor are now done. Next items of business are:
1. get raw and accepted vertices to Megan. She will begin with the raw 1.0 files, import as a pandas dataframe. Then she'll combine the "raw 2.0" vertices (raw, but with primary [and secondary] vertices corrected). THEN, she'll correlate raw vertices with accepted vertices.
2. Would be good to find a way for DSelector to output ONLY the TRUTH track permutation so that we can check vertexing. Not sure if this would be a good approach for training NNs, too.
UPDATE: There is a Get_IsTrueCombo() member of ParticleCombo, but it doesn't work for displaced vertices (see google group). Instead, it looks like I can check each charged particle hypothesis' thrown particle index. What I should do is add another branch to the output tree that is a boolean of whether the thrown topo is matched.
3. Need to understand how vertexing works. Does the kinematic fitter try to match the p mu- vertex with the Lambda flight path from the primary vertex? If so, this might NOT be helpful, as it could make fastpi events look too much like sl_mu events. So, can we access the pre-KF Lambda vertex in DSelector? Vertexing happens at hd_root time, right? Call Simon?
4. Would be good to include detector information (timing, dEdx, etc) in the data stream that goes into NNs. Find out how to do this in DSelector.
5. Look for conferences next year, and get Megan and Kripa on track for a poster presentation.
1. get raw and accepted vertices to Megan. She will begin with the raw 1.0 files, import as a pandas dataframe. Then she'll combine the "raw 2.0" vertices (raw, but with primary [and secondary] vertices corrected). THEN, she'll correlate raw vertices with accepted vertices.
2. Would be good to find a way for DSelector to output ONLY the TRUTH track permutation so that we can check vertexing. Not sure if this would be a good approach for training NNs, too.
UPDATE: There is a Get_IsTrueCombo() member of ParticleCombo, but it doesn't work for displaced vertices (see google group). Instead, it looks like I can check each charged particle hypothesis' thrown particle index. What I should do is add another branch to the output tree that is a boolean of whether the thrown topo is matched.
3. Need to understand how vertexing works. Does the kinematic fitter try to match the p mu- vertex with the Lambda flight path from the primary vertex? If so, this might NOT be helpful, as it could make fastpi events look too much like sl_mu events. So, can we access the pre-KF Lambda vertex in DSelector? Vertexing happens at hd_root time, right? Call Simon?
4. Would be good to include detector information (timing, dEdx, etc) in the data stream that goes into NNs. Find out how to do this in DSelector.
5. Look for conferences next year, and get Megan and Kripa on track for a poster presentation.
Tuesday, May 21, 2019
MC is now working, back to DSelector
The changes that I made to the ccdb call have paid off. Of the 900 MC files that I submitted earler, all but 3 have completed within three hours, and only about 12-15 of those that completed exited without finishing. Resubmitting these.
Next up: DSelector and vertices for Megan.
Recall how to run a DSelector over a TChain of root trees:
ifarm1402.jlab.org>
ifarm1402.jlab.org> root
root [0] TChain chain("sldmuon_p4vert_Tree");
root [1] chain.Add("../mc_files/kL_acc_mc/sl_mu/root/sl_mu_tree_sldmuon_p4vert_1*")
(Int_t) 293
root [2] .x $ROOT_ANALYSIS_HOME/scripts/Load_DSelector.C
root [3] sldmuon_p4vert_Tree->Process("DSelector_sld_muon_p4only.C+")
Next up: DSelector and vertices for Megan.
Recall how to run a DSelector over a TChain of root trees:
ifarm1402.jlab.org>
ifarm1402.jlab.org> root
root [0] TChain chain("sldmuon_p4vert_Tree");
root [1] chain.Add("../mc_files/kL_acc_mc/sl_mu/root/sl_mu_tree_sldmuon_p4vert_1*")
(Int_t) 293
root [2] .x $ROOT_ANALYSIS_HOME/scripts/Load_DSelector.C
root [3] sldmuon_p4vert_Tree->Process("DSelector_sld_muon_p4only.C+")
mc seems to be running correctly now... [NOPE]
Changed to -PJANA:MAX_RELAUNCH_THREADS=20 (i.e., from 10 to 20) in mcsmear call. Ten files running now with minimal difference between wall and cpu times.
Update: Nope, they failed... Still failing at the mcsmear stage.
On the cluster: 10/10 sl_mu jobs (only 1000 events each) crashed.
Checking output files. One difference is the valueof JANA_CALIB_URL.
On ifarm --> mysql://ccdb_user@hallddb.jlab.org/ccdb [works...]
On cluster --> sqlite:////work/halld/ccdb_sqlite/9/ccdb.sqlite [doesn't work???]
Adding to the sl_mu auger file:
setenv JANA_CALIB_URL mysql://ccdb_user@hallddb.jlab.org/ccdb
Rerunning 10 files.
OK, this seems to be the problem. I have no idea why it is intermittent, but I gave up on learning mysql looooong ago (2006).
Upping to 10k events and rerunning all mc.
Update: Nope, they failed... Still failing at the mcsmear stage.
On the cluster: 10/10 sl_mu jobs (only 1000 events each) crashed.
Checking output files. One difference is the valueof JANA_CALIB_URL.
On ifarm --> mysql://ccdb_user@hallddb.jlab.org/ccdb [works...]
On cluster --> sqlite:////work/halld/ccdb_sqlite/9/ccdb.sqlite [doesn't work???]
Adding to the sl_mu auger file:
setenv JANA_CALIB_URL mysql://ccdb_user@hallddb.jlab.org/ccdb
Rerunning 10 files.
OK, this seems to be the problem. I have no idea why it is intermittent, but I gave up on learning mysql looooong ago (2006).
Upping to 10k events and rerunning all mc.
Monday, May 20, 2019
mc toolchain working on nodes, not cluster
I have the mc toolchain working on the farm nodes. Processed a few batches of 100 sl_mu events with p4vert, p4only and p pi- reconstruction hypotheses. The root trees output look good.
I also rebuilt my local halld_recon, which took several hours.
Goal for today is to understand why mc recently quit working on the cluster. Recent output files show that hdgeant4 seems to be working OK (produces a reasonably sized output file). Next step in chain is mcsmear, and then hd_root.
UPDATE: I think that I fixed most of the problem. Added flags to mcsmear call in kL_mc_job_G4.py (and similar for ppim_fastpi script):
smear_cmd = 'mcsmear -PJANA_CALIB_CONTEXT="variation=mc" -PJANA:MAX_RELAUNCH_THREADS=10 '
Code now runs for sl_mu on the cluster. Deleting all previous MC files and rerunning.
I also rebuilt my local halld_recon, which took several hours.
Goal for today is to understand why mc recently quit working on the cluster. Recent output files show that hdgeant4 seems to be working OK (produces a reasonably sized output file). Next step in chain is mcsmear, and then hd_root.
UPDATE: I think that I fixed most of the problem. Added flags to mcsmear call in kL_mc_job_G4.py (and similar for ppim_fastpi script):
smear_cmd = 'mcsmear -PJANA_CALIB_CONTEXT="variation=mc" -PJANA:MAX_RELAUNCH_THREADS=10 '
Code now runs for sl_mu on the cluster. Deleting all previous MC files and rerunning.
Saturday, May 18, 2019
ppim running, sl_mu not...
MC update: ppim MC is (mostly) completing on the cluster, though wall time is still much larger than CPU time. Not sure what the issue is, but am nearly done with ppim MC. Will take a look at it this week.
sl_mu MC is NOT running. err files look like this
sl_mu MC is NOT running. err files look like this
Looking at the output file, looks like hdgeant4 is running and exiting OK, but mcsmear is writing 0 events to the hddm file. Going to submit a few more MC runs with an extra ls in them to sort this out.ls -lWarning from GlueXDetectorConstruction::ConstructSDandField - unsupported sensitive volume TAC1 found in geometry definiti\on.G4WT0 > Warning from GlueXDetectorConstruction::ConstructSDandField - unsupported sensitive volume TAC1 found in geometry defin\ition.JANA ERROR>>thread 0 has stalled on run:31000 event:101900001JANA ERROR>> Thread 0 (7f4082a8c700) hasn't responded in 31 seconds. (run:event=31000:101900001) Cancelling ...JANA ERROR>> Caught HUP signal for thread 0x7f4082a8c700 thread exiting...JANA ERROR>>JANA ERROR>> Last thread to lock output file mutex: 0x7f4082a8c700JANA ERROR>> Attempting to unlock mutex to avoid deadlock.JANA ERROR>> However, the output file is likely corrupt atJANA ERROR>> this point and the process should be restarted ...JANA ERROR>>JANA ERROR>>JANA ERROR>> Automatic relaunching of threads is disabled. If you wish toJANA ERROR>> have the program relaunch a replacement thread when a stalledJANA ERROR>> one is killed, set the JANA:MAX_RELAUNCH_THREADS configurationJANA ERROR>> parameter to a value greater than zero. E.g.:JANA ERROR>>JANA ERROR>> jana -PJANA:MAX_RELAUNCH_THREADS=10JANA ERROR>>JANA ERROR>> The program will quit now.cat: /usr/local/etc/motd: No such file or directorystty: standard input: Inappropriate ioctl for device/var/spool/torque-4.2.7/mom_priv/jobs/51868478.farmpbs14.SC: line 49: 999 Terminated nice -n 5 /scratch/auger/$J\OB_ID/command.sh >> $WORK_DIR/job.out 2>> $WORK_DIR/job.err=>> PBS: job killed: walltime 86443 exceeded limit 86400
Wednesday, May 15, 2019
How to get measured (i.e. from tracking, pre-KF) vertex from hd_root output?
Poking around in output root trees while I wait for MC to run.
(Which, btw, seems like there are some disk slowdowns at the lab and lots of my MC is timing out. Hitting walltimes of 13+ hrs. for jobs that should only take ~2 hrs.)
I can find the pre- and post-KF K+/beam vertex. However, it looks like only the post-KF vertex is available for Decaying_Lambda, proton, and muon.
I confess that I don't understand the ChargedHypo branch. Will keep digging.
MAYBE:
Thrown_X4: contains x4s for all of the thrown particles (makes sense with numbers of tracks, 5600 "events" match with 1400 events that make it through recon)
ChargedHypo_X4_Measured: contains x4s for all measured charged particles. 5177 events. Distribution looks mostly like the target shape, diminished Lam decay tails... Not sure what's up with this. Is recon trying to vertex each charged particle with the beam, but then KF gets the vertexing combos right???
Need to ask when the "measured" vertexing is happening for p and mu-.
(Which, btw, seems like there are some disk slowdowns at the lab and lots of my MC is timing out. Hitting walltimes of 13+ hrs. for jobs that should only take ~2 hrs.)
I can find the pre- and post-KF K+/beam vertex. However, it looks like only the post-KF vertex is available for Decaying_Lambda, proton, and muon.
I confess that I don't understand the ChargedHypo branch. Will keep digging.
MAYBE:
Thrown_X4: contains x4s for all of the thrown particles (makes sense with numbers of tracks, 5600 "events" match with 1400 events that make it through recon)
ChargedHypo_X4_Measured: contains x4s for all measured charged particles. 5177 events. Distribution looks mostly like the target shape, diminished Lam decay tails... Not sure what's up with this. Is recon trying to vertex each charged particle with the beam, but then KF gets the vertexing combos right???
Need to ask when the "measured" vertexing is happening for p and mu-.
MC still runinng
Most (~80%) of the jobs from yesterday timed out. Resubmitted today and they seem to be running smoothly.
Tuesday, May 14, 2019
Monday, May 13, 2019
Getting back into analysis! Status update.
This blog will serve as a "lab notebook" for my analysis of hyperon semi-leptonic decays using the GlueX experiment.
Today is the end of the Spring 2019 semester, and I'm getting back into the analysis. Current issues are with vertices:
1. A non-trivial number of signal events have both primary and lambda vertices assigned to the same point in space (roughly 800 events out of 250k that make it through the analysis chain). This should not be possible. One explanation could be that the kinfit does not converge for these events (and I'm letting them into the datastream anyway???), and it freaks out and assigns all particles the same vertex location. This is speculation.
2. The vertexing for the Lambda decay vertex seems to be too good. There is no appreciable difference between the transverse discrepancy distribution for signal and fastpi background.
WHAT AM I DOING ABOUT IT NOW?
1. rebuilding fastpion hdgeant4. I think that Richard has updated this over the last few months, and a new build is good. On this point, I have a script that adds all of the necessary changes to the hdgeant4 that is gotten from git, but the script does NOT softlink the patches. Remember to
mmcc @ HDGeant4$ ln -s G4.10.02.p02fixes G4fixes
2. Reprocess all MC. Might want to check with Hao about the KF fix -- should I go back to using the main branch rather than his fork?
3. Reinvestigate how DSelector gets the vertices for particles... also what cuts I'm applying. Probably want to keep them as loose as possible at this stage.
4. Look at distributions of raw and processed vertices. Megan Perry will assist later in the summer.
K THX BAI
- - - -
Update -- probably want to pull the master HDGeant4 and modify it. Looks like a lot has been added, and the changes have not trickled into the primarypionzapper branch.
Update -- successfully rebuilt master... testing now... Had to change run number in control.in to 31000 (or presumably something else).
../bin/.../hdgeant4 run.max works, but there is a huge spew (~10k lines) of particle def output at the top of the output. Possibly find a way to disable this...
Now adding PPZ stuff.
Update -- wrapped decay now works, but (oddly) needed to comment out the BEAM card in control.in in order to get the particle gun card to work. I don't remember having to do this before...
hdgeant4 now give about 10k lines of G4ParticleDefinition output when it's called. Did some googling, but it seems that there is not an easy way to silence this. Someone from Daya Bay has a very detailed workaround for their toolchain. Will ask around the collaboration.
- - - -
Purging old MC and regenerating now. 300 * 10k = 3M events for each decay.
Today is the end of the Spring 2019 semester, and I'm getting back into the analysis. Current issues are with vertices:
1. A non-trivial number of signal events have both primary and lambda vertices assigned to the same point in space (roughly 800 events out of 250k that make it through the analysis chain). This should not be possible. One explanation could be that the kinfit does not converge for these events (and I'm letting them into the datastream anyway???), and it freaks out and assigns all particles the same vertex location. This is speculation.
2. The vertexing for the Lambda decay vertex seems to be too good. There is no appreciable difference between the transverse discrepancy distribution for signal and fastpi background.
WHAT AM I DOING ABOUT IT NOW?
1. rebuilding fastpion hdgeant4. I think that Richard has updated this over the last few months, and a new build is good. On this point, I have a script that adds all of the necessary changes to the hdgeant4 that is gotten from git, but the script does NOT softlink the patches. Remember to
mmcc @ HDGeant4$ ln -s G4.10.02.p02fixes G4fixes
2. Reprocess all MC. Might want to check with Hao about the KF fix -- should I go back to using the main branch rather than his fork?
3. Reinvestigate how DSelector gets the vertices for particles... also what cuts I'm applying. Probably want to keep them as loose as possible at this stage.
4. Look at distributions of raw and processed vertices. Megan Perry will assist later in the summer.
K THX BAI
- - - -
Update -- probably want to pull the master HDGeant4 and modify it. Looks like a lot has been added, and the changes have not trickled into the primarypionzapper branch.
Update -- successfully rebuilt master... testing now... Had to change run number in control.in to 31000 (or presumably something else).
../bin/.../hdgeant4 run.max works, but there is a huge spew (~10k lines) of particle def output at the top of the output. Possibly find a way to disable this...
Now adding PPZ stuff.
Update -- wrapped decay now works, but (oddly) needed to comment out the BEAM card in control.in in order to get the particle gun card to work. I don't remember having to do this before...
hdgeant4 now give about 10k lines of G4ParticleDefinition output when it's called. Did some googling, but it seems that there is not an easy way to silence this. Someone from Daya Bay has a very detailed workaround for their toolchain. Will ask around the collaboration.
- - - -
Purging old MC and regenerating now. 300 * 10k = 3M events for each decay.
Subscribe to:
Comments (Atom)