mixer.cc

Go to the documentation of this file.

#include "mixer.hh"
#include <iostream>

bool compare_eventhits (EventHit* first, EventHit* second)
{
  //A check for null would be appropriate here, but what would you do when you found null?
  return first->time<=second->time;
}

void verify_input(int argc, char **argv)
{
   if (argc < 6 || argc % 2 == 1)
    {
      fprintf (stderr,"Usage: %s timeRandmax nTubes <outfile.txt> n*[<infile_n.txt> probability_to_selectfrom_this] \n",argv[0]);
      exit(1);
    }
}

void setup_eventsToMix(int argc, char **argv)
{
  double cumSum = 0; //For norming the probabilities.
  for(int i = 2; i<(argc/2); i++)
    {
      eventsToMix.push_back(make_pair(new FileEventStream(argv[2*i], NumberOfTubes),atof(argv[2*i+1])));
      cumSum+=eventsToMix.back().second;
    }
  //now norm the probabilities.
  for(vector<pair<FileEventStream*, double> >::iterator it = eventsToMix.begin(); it!=eventsToMix.end(); it++)
    {
      it->second/=cumSum;
    }
}

FILE * setup_outputFile(char* filename)
{
  FILE * outFile = fopen(filename, "w");
  if(outFile==NULL)
    {
      fprintf (stderr,"Could not open output file %s, exiting.", filename);
      exit(1);
    }
  return outFile;
}

void print_to_outfile(EventHit * it, FILE * outFile)
{
  fprintf(outFile,"%d %d %d %.15e ", it->eventno, it->Z, it->A, it->time);
 
  if(it->UFSP < DOUBLE_EPS)
    fprintf(outFile,"0 ");
  else
    fprintf(outFile,"%.15e ",it->UFSP);
  
  if(it->LFSP<DOUBLE_EPS)
    fprintf(outFile,"0 ");
  else
    fprintf(outFile,"%.15e ",it->LFSP);
  
  if(it->BSP<DOUBLE_EPS)
    fprintf(outFile,"0 ");
  else
    fprintf(outFile,"%.15e ",it->BSP);
  
  if(it->Ge<DOUBLE_EPS)
    fprintf(outFile,"0");
  else
    fprintf(outFile,"%.15e",it->Ge);
  
  for(int j = 0; j<NumberOfTubes; j++)
    {
      if(it->DET[j]<DOUBLE_EPS)
        fprintf(outFile," 0");
      else
        fprintf(outFile," %.15e", it->DET[j]);
    }
  fprintf(outFile,"\n");
}


void do_mixing(double randMax, FILE * outFile)
{
  double currentTime = 0;
  list<EventHit*> buffer;
  bool endLoop = false;
  long long nPrinted = 0;
  do
    {
      //Pick a new event with weighted probabilities for the different files.
      currentTime+=randMax*((double)rand()/(double)RAND_MAX);
      double u01 = ((double)rand()/(double)RAND_MAX);
      for(vector<pair<FileEventStream*, double> >::iterator it = eventsToMix.begin(); it!=eventsToMix.end(); it++)
        {
          u01-=it->second;
          if(u01<=0)
            {
              if(it->first->hasMoreEvents())
                {
                  Event * tempEvent = it->first->getNextEvent();
                  for(vector<EventHit*>::iterator it = tempEvent->hit.begin(); it!=tempEvent->hit.end(); it++)
                    {
                      (*it)->time+=currentTime;
                      buffer.push_back(*it);
                    }
                }
              else
                {
                  endLoop=true;
                }
              break;
            }
        }

      list<EventHit*> toPrint;

      for(list<EventHit*>::iterator it = buffer.begin(); it!=buffer.end(); it++)
        {
          if((*it)->time<currentTime || endLoop)
            {
              toPrint.push_back(*it);
              it = buffer.erase(it);
              --it;
            }
        }
      
      toPrint.sort(compare_eventhits);

      for(list<EventHit*>::iterator it = toPrint.begin(); it!=toPrint.end(); it++)
        {
          print_to_outfile(*it, outFile);
          ++nPrinted;
          delete *it;
        }

      toPrint.clear();
      

    }while(!endLoop);
}



int main(int argc, char **argv)
{
  verify_input(argc, argv);
  srand ( time(NULL) );
  
  double randMax = atof(argv[1]);
  NumberOfTubes = atoi(argv[2]);
  setup_eventsToMix(argc, argv);
  FILE * outFile = setup_outputFile(argv[3]);

  //Everything should be set up for the actual mixing now.
  do_mixing(randMax, outFile);
  
  //Clean up.
  fclose(outFile);
  return 0;
}