source: trunk/MagicSoft/Simulation/Corsika/MaxiSingle/maxisingle.cxx@ 14809

/////////////////////////////////////////////////////////////////////////////
//
//  PROGRAM MaxiSingle
//
//  MAXI files (multiple event) built form SINGLE event files
//  version 1.0
//
/////////////////////////////////////////////////////////////////////////////
//
//
//  AUTHOR Carles Domingo
//  May 15, 2001
//
/////////////////////////////////////////////////////////////////////////////
//
//  PROCEDURE
//
//  This program reads a set of files in the old MMCS output format (one
//  file per event) and writes the information from those files in the
//  new MMCS output format (one file per run).
//
//  To perform this task, the information to fill in the Run headers, the
//  run end records and the event end record (which were missing in the
//  old format files) has to be rebuilt from the information read for
//  the events.
//
//  Some new CORtype classes had to be defined to cope with this 
//  information in an adequate way. The classes COREventEnd, CORRunHeader
//  and CORRunEnd are now defined and operational.
//
//  In order to run the program you must supply as an argument the full
//  path of a directory containing the output files (cer and dat), starting
//  from 000001 of a given (old) MMCS run. The program will automatically  
//  create inside this directory two files (named cer_allfiles and 
//  dat_allfiles) containing the information of the complete run in the
//  new single file per run output format. Afterwards you can use this files
//  for the new reflector, camera, and guisimone programs. 
//
/////////////////////////////////////////////////////////////////////////////
//
//  REMARKS
//
//  Tricky stuff for locating the event ends in the dat FILES!! See comments
//  on the program
//
//
//   ATTENTION: The estructure of the input CER files and that of 
//   the DAT files is different in the "one file per event" version 
//   of the MMCS program. 
//
//   CER files have:
//     Event header EVTH
//     A bunch of cerenkov photons (in blocks of 7*39=273 bytes)
//
//   while DAT files have:
//     Event header EVTH
//     A bunch of particles (in blocks of 7*39=273 bytes)
//     Event end EVTE
//
/////////////////////////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////// 
//
//  NOTE
//
//  If anybody modifies this program, make sure it is well documented !!!!
//
/////////////////////////////////////////////////////////////////////////////


#include <stdlib.h>
#include <stdio.h>

#include "COREventHeader.hxx"
#include "COREventEnd.hxx"
#include "CORParticle.hxx"
#include "CORRunHeader.hxx"
#include "CORRunEnd.hxx"

int main(int argc, char **argv)
{
  char  incername[120] ;
  char  indatname[120] ;
  char  outcername[120] ;
  char  outdatname[120] ;

  Float_t  NumberRun ;
  Float_t  EvtsProc ;  
  Float_t  TestVal ;

  ifstream   incerfile ;
  ifstream   indatfile ;
  ofstream   outcerfile ;
  ofstream   outdatfile ; 
  
  COREventHeader Event ;
  COREventHeader EventDat ;
  COREventEnd    EventE ;   
  CORParticle    Photon ;
  CORParticle    Particle ; 
  CORRunHeader   Run ;
  CORRunEnd      RunE ;
  
  cout << " ============================================================" << endl ; 
  cout << " MaxiSingle" << endl ; 
  cout << "       " << endl ; 
  cout << " MAXI (multiple events) file from SINGLE event multiple files" << endl ; 
  cout << "       " << endl ; 
  cout << " ============================================================" << endl ; 

  if ( argc <= 1 ) {
    cout << endl ; 
    cout << " INFO: You have to start the program with "<< endl << endl ; 
    cout << "        -> maxisingle DIRECTORY_WITH_CER_FILES" << endl << endl ; 
    cout << "        no SLASH at the end of the directory name"<< endl ; 
    cout << "        (example: -> maxisingle /hd123/Protons   "<< endl ;
    cout << "  " << endl ; 
    exit (-1) ; 
  }

  //   create and open the OUTput filenames
  
  sprintf ( outcername, "%s/cer_allfiles", argv[1] ) ;
  sprintf ( outdatname, "%s/dat_allfiles", argv[1] ) ;
  outcerfile.open (outcername) ;
  if ( outcerfile.bad() ) { 
    cout << "Cannot open CER output file:  " << outcername << endl ;
    exit (1) ; 
  } 
  outdatfile.open (outdatname) ;
  if ( outdatfile.bad() ) { 
    cout << "Cannot open DAT output file:  " << outdatname << endl ;
    exit (1) ; 
  }   

  EvtsProc = 0 ;
  for (int i_cer = 1; i_cer <= 100000; i_cer++ ) {

    //   inform about progress
      
    if (!( i_cer %100) ) cout << " Processing event number: " << i_cer << endl ; 
   
    //   create the file names 
      
    sprintf ( incername, "%s/cer%06d", argv[1], i_cer ) ;
    sprintf ( indatname, "%s/dat%06d", argv[1], i_cer ) ; 
      
    //   try to open the files
      
    incerfile.open( incername );
    if ( incerfile.bad() )  break ;
    indatfile.open( indatname );
    if ( indatfile.bad() )  break ;

    //   read and write the event header
      
    Event.read( incerfile ) ;
    EventDat.read( indatfile ) ;
  
    //   if first event, create and write a RUN header in
    //   both outcerfile and outdatfile
  
    if (i_cer == 1 ) { 
      //      Event.print() ;
      Run.transport (&Event) ;
      NumberRun = Run.get_number () ;
      Run.write (outcerfile) ;
      Run.write (outdatfile) ;
      //      Run.print () ;
    }
   
    Event.write( outcerfile ) ;
    EventDat.write( outdatfile ) ;
    EvtsProc++ ;
      
    //   loop over the photons in incerfile, read and write them
      
    while( ! (incerfile.eof() || incerfile.bad() )) {
      Photon.read ( incerfile ) ;
      Photon.write ( outcerfile ) ;
    }
    
    //   loop over the particles in indatfile, read and write them
    //   This procedure is somewhat tricky, in order to find the event end record
    
    while( ! (indatfile.eof() || indatfile.bad() )) {
      
      Particle.read ( indatfile ) ; 
      //      Particle.print () ;       

      //  TRICKY STUFF
      //
      //  here we get the value which is in the first position of the particle
      //  record. Whenever there is a particle, this Float_t variable has an
      //  integer value. If there is no particle (but we still have a particle 
      //  type record that MMCS has filled with zeroes), the value is 0. But if
      //  the EVENT END record is found, we find the Char_t identifier 'EVTE' 
      //  in that first position. If 'EVTE' is read as Float_t it returns a value 
      //  of 3397.39 both in OSF1 and Linux, so we can easyly check for that value,
      //  also making sure that it is not integer (any number between 3397.38 and
      //  3397.40 CANNOT be integer)... so it is sure it is NOT a real particle.

      TestVal = Particle.get_id () ;
      if (TestVal > 3397.38 && TestVal < 3397.40) {

        //  Go back in the file to read again the record not as a particle record
        //  but as EVENT END

        indatfile.seekg( -28, ios::cur  ) ;     
        EventE.read( indatfile );
        break ; 
      }  
      Particle.write ( outdatfile ) ;
      //      Particle.print () ;
    }

    //   close input files and write the EVENT end record in outputs
    
    incerfile.close () ;
    indatfile.close () ;
    // EventE.print () ;
    EventE.write( outcerfile );
    EventE.write( outdatfile );
    
  } 
    
  //  write the RUN end record and close output files      
        
  cout << " ---> Processed a total of " << EvtsProc << " events" << endl ;
  RunE.fill ( NumberRun , EvtsProc ) ;
  RunE.write( outcerfile );
  RunE.write( outdatfile );
  outcerfile.close ();
  outdatfile.close ();
  
}