fillSpaceChargeMaps/macros/simulateTimestamps.C

0001 #include <TTree.h>
0002 #include <TFile.h>
0003 #include <gsl/gsl_randist.h>
0004 #include <gsl/gsl_rng.h>
0005
0006 #include <stdio.h>
0007 #include <stdlib.h>
0008
0009 // cppcheck-suppress unknownMacro
0010 R__LOAD_LIBRARY(libgslcblas.so)
0011 R__LOAD_LIBRARY(libgsl.so)
0012
0013 //_______________________________________________________________
0014 void simulateTimestamps()
0015 {
0016
0017   // random generator
0018   auto rng = gsl_rng_alloc(gsl_rng_mt19937);
0019   //gsl_rng_set( rng, time(0) );//time(0) is the seed
0020   gsl_rng_set( rng, 1 );
0021
0022   // time interval (s) between two bunch crossing
0023   /* value copied from generators/phhepmc/Fun4AllHepMCPileupInputManager.cc */
0024   static constexpr double deltat_crossing = 106e-9;
0025
0026   // triggered rate (Hz)
0027   // static constexpr double trigger_rate = 1e5; // Hijing 100kHz
0028   static constexpr double trigger_rate = 5e4; // Hijing 50kHz
0029   //static constexpr double trigger_rate = 3e6; // pp collisions
0030
0031   // mean number of collision per crossing
0032   static constexpr double mu = trigger_rate*deltat_crossing;
0033
0034   // print configuration
0035   printf("SimulateCollisions - deltat_crossing: %f \n", deltat_crossing );
0036   printf("SimulateCollisions - trigger_rate: %f \n" , trigger_rate);
0037   printf("SimulateCollisions - mu: %f \n" , mu);
0038
0039   // number of requested triggers (should correspond approximately to 1/10 second of data taking)
0040   // static constexpr uint ntrigtot = 1e4;
0041   //static constexpr uint ntrigtot = 6e6;
0042   static constexpr uint ntrigtot = 1e6;
0043   // write all timestamps to a file
0044   FILE *fptr;
0045   fptr = fopen("./data/timestamps_50kHz_1M.txt","w");
0046   fprintf(fptr,"// bunchcrossin id; time (ns)\n");
0047   fprintf(fptr,"// assuming 106ns between bunches\n" );
0048
0049
0050   // running collision time
0051   int64_t bunchcrossing = 0;
0052   float time = 0;
0053   float timeD = 0;
0054   TTree *timestamps = new TTree("timestamps","beamcrossings timestamps");
0055   timestamps -> Branch("bc",  &bunchcrossing,   "bc/I");
0056   timestamps -> Branch("t",  &time,   "t/F");
0057   timestamps -> Branch("dt",  &timeD,   "dt/F");
0058
0059   // keep track of the last collision time
0060   double previous_trigger_time = 0;
0061
0062   // generate triggers
0063   for( int itrig = 0; itrig < ntrigtot; )
0064   {
0065     ++ bunchcrossing;
0066     time += deltat_crossing;
0067
0068     auto ntrig = gsl_ran_poisson( rng, mu );
0069     for( uint i = 0; i < ntrig; ++i )
0070     {
0071       fprintf(fptr,"%d %d\n",int(bunchcrossing), int(time*1e9) );
0072       timeD = (time-previous_trigger_time)*1e9;
0073       timestamps ->Fill();
0074       previous_trigger_time = time;
0075       ++itrig;
0076       if( itrig%100000==0 ) printf("SimulateCollisions - itrig: %d \n", itrig );
0077     }
0078   }
0079
0080   // printout last trigger time
0081   printf("SimulateCollisions - last trigger time: %f \n", time );
0082
0083   fclose(fptr);
0084   //out.close();
0085   gsl_rng_free(rng);
0086
0087   TFile *outf = new TFile("./data/timestamps_50kHz_1M_t.root","recreate");
0088   timestamps -> Write();
0089
0090   outf -> Write();
0091   outf -> Close();
0092
0093 }