Take a look at: https://www.cs.rice.edu/~vo9/deeplearning/2019/slides/Rivanna.pdf
Here are the updated instruction of generating Monte Carlo for SpinQuest simulations and analysis
- Login to Rivanna (https://ood.hpc.virginia.edu/pun/sys/dashboard)
- Create a folder to submit your jobs from on /project/ptgroup/ or /project/UVA-Spin
- The scripts to generate Monte Carlo (MC) events on Rivanna are in:
/project/ptgroup/work/MC_Generation
You can make more and store the generating scripts here but label them in a way that is easy to understand. The present naming convention is <Exp>_<channel>_<vertex origin>- Option for <Exp> are 906 or 1039. The Geant4 geometry is different for the setup of E906 and E1039, so for SpinQuest make sure you are using the E1039 geometry.
- Option for <channel> are : Drell-Yan (DY), Jpsi (JPsi), Psi' (Psip) Pion background (Pion), Kaon background (Kaon), Random-multi muons (MultiMuon), Single Muons + or- (SingMup/m)
- Option for <Vertex origin> : target (Target), Beam Dump (Dump), Everything that seen by the beam (All), The gap between the target and dump (TargetDumpGap), and arbitrary vertex origin (Manual, should specify numbers such as x,y,z→25,25,100 would be a target size with diameter 25cm and 100cm long.
- You can copy any MC scripts to your directory to submit. The scripts to generate Monte Carlo (MC) events on Rivanna are in:
/project/ptgroup/work/MC_Generation
Example of copying: cd /project/ptgroup/YourFolder
rsync -av --exclude scratch /project/ptgroup/work/MC_Generation/DY_Target_script .
We will keep track of all MC generation scripts so we can all see how a particular MC batch was produced. Details of each one should be kept in ReadMe in that directory. - Other MC configurations you make that are not included in this directory you should copy over and make a note about so other people can use it as well.
- To summit you jobs to produce the MC navigate to the directory where you want to run the script.
- Setup the environment using: source /project/ptgroup/spinquest/this-e1039.sh command
Submit your job using the following command.
./jobscript.sh <Outputfoldername> <Number of jobs> <Number of Events per Job> For example: ./jobscript.sh DY_Target_1M 100 10000 (Here we need to generate ~1M events. The accepted events depend on the channel). Strongly recommended not to exceed 10K events per job.
The jobscript.sh should look like this:
#!/bin/bash
dir_macros=$(dirname $(readlink -f $BASH_SOURCE))
jobname=$1
njobs=$2
nevents=$3
echo "njobs=$njobs"
echo "nevents=$nevents"
work=/scratch/$USER/MC/$jobname
mkdir -p $work
chmod -R 01755 $work
cd $dir_macros
for (( id=1; id<=$njobs; id++ ))
do
mkdir -p $work/$id/
chmod -R 01755 $work/$id
cd $work/$id/
cp $dir_macros/*.C .
cp $dir_macros/*.cfg .
cp $dir_macros/*.txt .
cp $dir_macros/*.slurm .
sed -i "s/1234/$nevents/" Fun4Sim.C
echo "submitting job number = $id"
sbatch grid.slurm
Which will submit the job according to the grid.slurm script:
#!/bin/sh
#SBATCH --ntasks=1
#SBATCH --ntasks-per-node=1
#SBATCH --time=6:00:00
#SBATCH --output=slurm.out
#SBATCH --error=slurm.err
#SBATCH --partition=standard
#SBATCH -A spinquest
export DISPLAY=:0.0
root -l Fun4Sim.C
root -l filter.C
- This is the slurm script for Rivanna (https://www.rc.virginia.edu/userinfo/hpc/slurm/)
- sbatch -A gives the Rivanna allocation. Both spinquest and spinquest_standard should be available.
- You can use "squeue -u user", to check the status of your jobs (or use the "Active Jobs" tab on your UVA OpenOnDemand web page: https://ood.hpc.virginia.edu/pun/sys/dashboard/activejobs?jobcluster=all&jobfilter=user).
- Navigate to "/project/ptgroup/script" and then implement the following command (this script by-default will skip bad-files/corrupted-files/etc.)
$ ./merge_mc_prod.sh /scratch/<your_MC_output_file_location>
This will create a folder inside /project/ptgrpup/spinquest/MC_merge_files" with the same name as your MC_output_folder - Copy the "merged_trackQA.root or merged_trackQA_v2.root" file from the folder that you've newly generated by ./merge_mc_prod.sh and copy it to "/project/ptgroup/spinquest/MC_storage" location. Rename the file into <Channel>_<Vertex Origin>_<Number of accepted event>.root. for example: JPsi_Dump_300K.root
- The output files of the submitted jobs can be found at: /sfs/weka/scratch/<username>/MC
Additional Information about changing the Setup
The script that you need to change for any particular configuration is the Fun4Sim.C
Most of this configuration file you will not need to change, but depending on what you want to do here are some pointers
Choice of the Generators:
This is how you would select what physics generator you want to use (for normal DY or J/Psi use pythia):
const bool gen_pythia8 = true;
const bool gen_cosmic = false;
const bool gen_gun = false;
const bool gen_particle = false;
const bool read_hepmc = false;
const bool gen_e906legacy = false;
gen_pythia8:
1. If you want to use Pythia8, you can load the config file using this function. pythia8->set_config_file("phpythia8_DY.cfg") or pythia8->set_config_file("phpythia8_Jpsi.cfg");
depending on your work.
2. In the configuration file, you can change the configuration based on your need. For example, you can set the mass cut at Pythia8 by putting this line: PhaseSpace:mHatMin = 4.0.
3. In addition to the pythia configuration file, you can add more variables using PHPy8ParticleTrigger. A basic setup has shwn below, but you can find more function here.
// pythia8
if(gen_pythia8) {
PHPythia8 *pythia8 = new PHPythia8();
//pythia8->Verbosity(99);
pythia8->set_config_file("phpythia8_DY.cfg");
if(SQ_vtx_gen) pythia8->enableLegacyVtxGen();
else{
pythia8->set_vertex_distribution_mean(0, 0, target_coil_pos_z, 0);
}
pythia8->set_embedding_id(1);
se->registerSubsystem(pythia8);
pythia8->set_trigger_AND();
PHPy8ParticleTrigger* trigger_mup = new PHPy8ParticleTrigger();
trigger_mup->AddParticles("-13");
//trigger_mup->SetPxHighLow(7, 0.5);
//trigger_mup->SetPyHighLow(6, -6);
trigger_mup->SetPzHighLow(120, 10);
pythia8->register_trigger(trigger_mup);
PHPy8ParticleTrigger* trigger_mum = new PHPy8ParticleTrigger();
trigger_mum->AddParticles("13");
//trigger_mum->SetPxHighLow(-0.5, 7);
//trigger_mum->SetPyHighLow(6, -6);
trigger_mum->SetPzHighLow(120, 10);
pythia8->register_trigger(trigger_mum);
}
}
gen_e906legacy:
This generator is adopted from the SeaQuest experiment. If you interested in generating dimuon mass, especially in the low mass regions, you can use the script, because generating
more statistics using gen_pythia8 is computationally much more expensive. And If you use gen_e906legacy generator, you will need to save the event weight information, since
each event has a generator weight. For example, the DY has been used below:
if(gen_e906dim){
SQPrimaryParticleGen *e906legacy = new SQPrimaryParticleGen();
const bool pythia_gen = false;
const bool drellyan_gen = true;
const bool JPsi_gen = false;
const bool Psip_gen = false;
if(drellyan_gen){
e906legacy->set_xfRange(0.1, 0.5); //[-1.,1.]
e906legacy->set_massRange(3.0, 10.0);
e906legacy->enableDrellYanGen();
}
if(Psip_gen){
e906legacy->set_xfRange(0.1, 0.5); //[-1.,1.]
e906legacy->enablePsipGen();
}
if(JPsi_gen){
e906legacy->set_xfRange(0.1, 0.5); //[-1.,1.]
e906legacy->enableJPsiGen();
}
}
Generating Single or Multi muons from using the multi particle gun
For making a single or multiple muons you can use multi particle gun. The portion of the script below can be found on the SimChainDev module. Set the gen_particle true if you would like to generate single muons from multi particle guns. Some user function are available on the doxygen page.
if(gen_particle) {
PHG4SimpleEventGenerator *genp = new PHG4SimpleEventGenerator("MUP");
//genp->set_seed(123);
genp->add_particles("mu+", nmu); // mu+,e+,proton,pi+,Upsilon
if (legacyVtxGen) genp->enableLegacyVtxGen();
else{
genp->set_vertex_distribution_function(PHG4SimpleEventGenerator::Uniform,
PHG4SimpleEventGenerator::Uniform,
PHG4SimpleEventGenerator::Uniform);
genp->set_vertex_distribution_mean(0.0, 0.0, -300);
genp->set_vertex_distribution_width(100.0, 100.0, 400.0);
genp->set_vertex_size_function(PHG4SimpleEventGenerator::Uniform); [For Gaussian distributions use "Gaus"]
genp->set_vertex_size_parameters(0.0, 0.0);
}
Basic Configurations for the simulations:
const double target_coil_pos_z = -300;
const int nmu = 1;
int embedding_opt = 0;
const bool legacy_rec_container = true;
const bool do_collimator = true;
const bool do_target = true;
const bool do_e1039_shielding = true;
const bool do_fmag = true;
const bool do_kmag = true;
const bool do_absorber = true;
const bool do_dphodo = true;
const bool do_station1DC = false; //station-1 drift chamber should be turned off by default
const double target_l = 7.9; //cm
const double target_z = (7.9-target_l)/2.; //cm
const int use_g4steps = 1;
const double FMAGSTR = -1.044;
const double KMAGSTR = -1.025;
Here you can change the muon to mu+ or mu-, or change the or change the vertex interaction volume, for example this one creates a cylinder 1m radius around the beamline from 700 cm upstream of the dump to 100 cm inside the dump.
- Beam profile
- The profile here means the distribution shape of beam protons in X and Y (= R).
- When "legacyVtxGen = true"
- The shape is "Gaussian" at "R < 5*sigma" and "1 / R" at "R >= 5*sigma" (cf. https://github.com/E1039-Collaboration/e1039-core/blob/master/generators/E906LegacyVtxGen/SQPrimaryVertexGen.C#L290 )
- The Gaussian width ("sigma" in cm) is defined by "SIGX_BEAM" and "SIGY_BEAM" in Fun4Sim.C.
- When "legacyVtxGen = false"
- The shape is defined by the following functions of event generators.
- "set_vertex_distribution_function()" sets the shape to "Uniform" or "Gaussian".
- "set_vertex_distribution_mean()" sets the mean of the distribution.
- "set_vertex_distribution_width()" sets the half width in case of "Uniform" or "sigma" in case of "Gaussian".
- Distribution of z-vertex
- The z-vertex means the position of each generated event in Z.
- When "legacyVtxGen = true"
- The distribution shape follows the probability of beam-material interaction, namely the product of beam intensity and material amount at Z, where the beam intensity is attenuated over Z.
When "VTX_GEN_MATERIAL_MODE" is set to "Target" (via
rc->set_CharFlag("VTX_GEN_MATERIAL_MODE", "Target")
), events are generated only on the target material.- "VTX_GEN_MATERIAL_MODE" can be set to "Dump" or "TargetDumpGap" to enable only the dump material or the material (i.e. air) in between the target and the dump.
- When "legacyVtxGen = false"
- The distribution shape is defined by "set_vertex_distribution_function()" etc. together with the beam profile.
- Target position
- The z-center and the length of the target is defined by "double target_coil_pos_z" and "double target_l" in Fun4Sim.C (-300 and -7.9 cm by default).
- "target_z" is the position of the target material relative to the target coils. It is 0 cm by default. If you change "target_l", you have to change "target_z" to "0" together, unless you intend to make the relative position non-zero.
- The x,y-center and the angle of the target are 0 by default. If you need to change them,
Copy "G4_Target.C" to your working directory;
cp /project/ptgroup/spinquest/core/default/macros/top/G4_Target.C .
- Edit your "G4_Target.C" file.
- Replace
#include <top/G4_Target.C>
with#include "G4_Target.C"
in "Fun4Sim.C".
To save the accepted tracks/Dimuons:
If you would like to save only accepted single muon or dimuon, you can use SQGeomAcc.
/// Save only events that are in the geometric acceptance.
SQGeomAcc* geom_acc = new SQGeomAcc();
geom_acc->SetMuonMode(SQGeomAcc::PAIR); // PAIR, PAIR_TBBT, SINGLE, SINGLE_T, etc.
geom_acc->SetPlaneMode(SQGeomAcc::HODO_CHAM); // HODO, CHAM or HODO_CHAM
geom_acc->SetNumOfH1EdgeElementsExcluded(4); // Exclude 4 elements at H1 edges
se->registerSubsystem(geom_acc);