/////////////////////////////////////////////////////////////////
//
// MCEventHeader_2
//
//  Created: Sat Oct 19 12:11:43 CEST 2002
//  Author:  Oscar Blanch Bigas
//  Purpose: Base class for reflector 0.6 EventHeader
//  Notes:   
//
/////////////////////////////////////////////////////////////////

// @T \newpage

// @section Source code of {\tt MCEventHeader_2.hxx}

/* @text
This section shows the include file {\tt MCEventHeader_2.hxx}
@endtext */

#ifndef MCEventHeader_2_Class
#define MCEventHeader_2_Class

// @subsection Include files

// @code
#include "TROOT.h"
#include "TObject.h"

#include <iostream.h>
#include <iomanip.h>
#include <fstream.h>
#include <stdlib.h>
#include <math.h>

#include "COREventHeader.hxx"

#define SIZE_OF_MCEVENTHEADER_2    181   /* floats */

// @endcode

// @subsection Class {\em MCEventHeader_2}: Definition

// @code
class MCEventHeader_2 {

protected:
  Float_t      EvtNumber;
  Float_t      PrimaryID;
  Float_t      Etotal;   
  Float_t      Thick0;   
  Float_t      FirstTarget;
  Float_t      zFirstInt;
  Float_t      p[3];     
  Float_t      Theta; 
  Float_t      Phi; 

  Float_t      NumRndSeq;
  Float_t      RndData[10][3];
  
  Float_t      RunNumber;
  Float_t      DateRun;
  Float_t      VersionPGM;

  Float_t      NumObsLev;
  Float_t      HeightLev[10]; 

  Float_t      SlopeSpec;
  Float_t      ELowLim;   
  Float_t      EUppLim;   

  Float_t	Ecutoffh;  
  Float_t	Ecutoffm;  
  Float_t	Ecutoffe;  
  Float_t	Ecutoffg;  
  
  Float_t	NFLAIN;
  Float_t	NFLDIF;
  Float_t	NFLPI0;
  Float_t	NFLPIF;
  Float_t	NFLCHE;
  Float_t	NFRAGM; 
  
  Float_t	Bx;
  Float_t	By;
  
  Float_t	EGS4yn;
  Float_t	NKGyn;
  Float_t	GHEISHAyn;
  Float_t	VENUSyn;
  Float_t	CERENKOVyn;
  Float_t	NEUTRINOyn;
  Float_t	HORIZONTyn;
  Float_t	COMPUTER;
  
  Float_t      ThetaMin;
  Float_t      ThetaMax;
  Float_t      PhiMin;
  Float_t      PhiMax;

  Float_t	CBunchSize;
  Float_t	CDetInX,CDetInY;
  Float_t	CSpacInX,CSpacInY;
  Float_t	CLenInX,CLenInY;
  Float_t	COutput;
  
  Float_t	AngleNorthX;
  Float_t	MuonInfo;
  
  Float_t	StepLength;
  Float_t	CWaveLower;       
  Float_t	CWaveUpper;       
  Float_t	Multipl;       
  Float_t	CorePos[2][20];   
  Float_t       SIBYLL[2];
  Float_t       QGSJET[2];
  Float_t       DPMJET[2];
  Float_t       VENUS_cross;
  Float_t       mu_mult_scat;
  Float_t       NKG_range;
  Float_t       EFRCTHN[2];
  Float_t       WMAX[2];
  Float_t       rthin_rmax;

  Float_t       viewcone_angles[2]; 
  /* (degrees) Inner and outer angles in 
   * Corsika's VIEWCONE option. This is
   * only possible with Corsika>6 versions. In
   * that case, PhiMin=PhiMax  and 
   * ThetaMin=ThetaMax (also in this header) 
   * indicate the axis of this cone.   
   */

  /* ^^^ Up to here, the info from the CORSIKA event header. */

  /* Telescope orientation: */
  Float_t	telescopePhi;    /* rad */
  Float_t	telescopeTheta;  /* rad */

  /* Time of first and last photon: */
  Float_t      TimeFirst;
  Float_t      TimeLast;

  /* 6 parameters and chi2 of the NKG fit to the longitudinal 
   * particle distribution (see CORSIKA manual for explanation): 
   */
  Float_t       longi_Nmax;
  Float_t       longi_t0;
  Float_t       longi_tmax;
  Float_t       longi_a;
  Float_t       longi_b;
  Float_t       longi_c;
  Float_t       longi_chi2;

  Float_t      CORSIKAPhs;   // Original photons written by Corsika
  Float_t      AtmAbsPhs;    // Photons absorved by the atmosphere
  Float_t      MirrAbsPhs;   // Photons absorved by the mirrors
  Float_t      OutOfMirrPhs; // Photons outside the mirror
  Float_t      BlackSpotPhs; // Photons lost in the "black spot"
  Float_t      OutOfChamPhs; // Photons outside the chamber
  Float_t      CPhotons;     // Photons reaching the chamber

    /* Now follow the fraction of photons reaching the camera produced by
     * electrons, muons and other particles respectively:
     */
  Float_t       elec_cph_fraction;
  Float_t       muon_cph_fraction;
  Float_t       other_cph_fraction;
  
  Float_t      dummy[7];     /* not used */

public:
  MCEventHeader_2(void) {}  // default constructor

  virtual ~MCEventHeader_2(void) {}  // default destructor
  
  // reads EventHeader from binary input stream
  Int_t read ( ifstream &is ) {
    is.read ( (char *)this, mysize() ); 
    return is.gcount();
  }

  // writes EventHeader to binary output stream
  Int_t write ( ofstream &os ) {
    os.write ( (char *)this, mysize() ); 
    return 0;
  }
  
  // get information about the EventHeader

  // get the event number
  inline Float_t get_evt_number (void) {return EvtNumber;}

  // get the primary type (GEANT code)
  inline Float_t get_primary ( void ) { return ( PrimaryID ); }

  // get the total primary energy
  inline Float_t get_energy ( void ) { return ( Etotal ); }

  // get Thick0
  inline Float_t get_thick0 (void) { return Thick0;}

  // get First Target
  inline Float_t get_first_target (void) {return FirstTarget;}

  // get z first interaction
  inline Float_t get_z_first_int (void) {return zFirstInt;}

  // get the initial zenith angle
  inline Float_t get_theta ( void ) { return ( Theta ); }

  // get the initial phi angle
  inline Float_t get_phi ( void ) { return ( Phi ); }

  // get the telescope zenith angle
  inline Float_t get_theta_CT ( void ) { return ( telescopeTheta ); }

  // get the telescope phi angle
  inline Float_t get_phi_CT ( void ) { return ( telescopePhi ); }

  // get the spectrum slope
  inline Float_t get_slope ( void ) { return ( SlopeSpec ); }

  // get height of first interaction (in cm)
  inline Float_t get_height ( void ) { return ( zFirstInt ); }

  // get the energy range of this run
  inline void get_energy_range ( Float_t *elow, Float_t *eup ) { 
    *elow = ELowLim;   
    *eup = EUppLim;
  }

  // get the theta range of this run
  inline void get_theta_range ( Float_t *thlow, Float_t *thup ) { 
    *thlow = ThetaMin;   
    *thup = ThetaMax;
  }

  inline void get_viewcone ( Float_t *rmin, Float_t *rmax){
    *rmin=viewcone_angles[0];
    *rmax=viewcone_angles[1];
}

  // get the energy range of this run
  inline void get_phi_range ( Float_t *plow, Float_t *pup ) { 
    *plow = PhiMin;   
    *pup = PhiMax;
  }

  // get the core position
  inline Float_t get_core ( Float_t *x, Float_t *y, Int_t ncore = 0 ) {
    *x = CorePos[0][ncore];
    *y = CorePos[1][ncore];
    return ( (Float_t) sqrt ( (*x)*(*x) + (*y)*(*y) ) );
  }

  // get the core position
  inline Float_t get_core ( Int_t ncore = 0 ) {
    return ( (Float_t) sqrt ( ((CorePos[0][ncore])*(CorePos[0][ncore]))+
                              ((CorePos[1][ncore])*(CorePos[1][ncore])) ) );
  }
  
  // get the core position in X 
  inline Float_t get_coreX ( Int_t ncore = 0 ) {
    return ( (Float_t) CorePos[0][ncore] );
  }

  // get the core position in Y 
  inline Float_t get_coreY ( Int_t ncore = 0 ) {
    return ( (Float_t) CorePos[1][ncore] );
  }

  //get photons at different stages
  inline Float_t get_CORSIKA (){ return CORSIKAPhs;}

  inline Float_t get_AtmAbs (){ return AtmAbsPhs;}

  inline Float_t get_MirrAbs (){ return MirrAbsPhs;}

  inline Float_t get_OutOfMirr (){ return OutOfMirrPhs;}

  inline Float_t get_BlackSpot (){ return BlackSpotPhs;}

  inline Float_t get_OutOfCham (){ return OutOfChamPhs;}

  inline Float_t get_CPhotons (){ return CPhotons;}

  inline Float_t get_NumObsLev (){return NumObsLev;}

  inline Float_t get_HeightLev (Int_t i) {return HeightLev[i];}

  inline Float_t get_VersionPGM () {return VersionPGM;}

  inline Float_t get_RunNumber () {return RunNumber;}

  inline Float_t get_DateRun () {return DateRun;}

  inline Float_t get_NumRndSeq () {return NumRndSeq;}

  inline Float_t get_CWaveLower() {return CWaveLower;}

  inline Float_t get_CWaveUpper() {return CWaveUpper;}

  inline Float_t get_ElecFraction() {return elec_cph_fraction;}

  inline Float_t get_MuonFraction() {return muon_cph_fraction;}

  inline Float_t get_OtherFraction() {return other_cph_fraction;}

  // transport from COREventHeader to MCEventHeader_2
  void transport ( COREventHeader *e );

  // write extreme times
  inline void put_times ( Float_t t1, Float_t t2 ) {
    TimeFirst = t1;
    TimeLast  = t2;
  }

  // get extreme times
  inline Float_t get_times ( Float_t *t1, Float_t *t2 ) {
    *t1 = TimeFirst;
    *t2 = TimeLast;
    return ( TimeLast - TimeFirst );
  }

  inline Float_t get_NKGfit(Float_t *Nmax, Float_t *t0, Float_t *tmax,
			    Float_t *a,    Float_t *b,  Float_t *c){
    *Nmax=longi_Nmax; *t0=longi_t0; *tmax=longi_tmax;
    *a=longi_a;       *b=longi_b;   *c=longi_c;
    
    return longi_chi2;

  }

  inline int mysize(void) 
    { return ( sizeof( float ) * SIZE_OF_MCEVENTHEADER_2 ); }

  // get deviations of the CT from the original shower direction
  inline Float_t get_deviations ( Float_t *t1, Float_t *t2 ) {
 
    float ct1,st1,cp1,sp1;
    float ct2,st2,cp2,sp2;
    float x1,y1,z1;
    float x2,y2,z2;

    ct1 = cos(Theta);
    st1 = sin(Theta);
    cp1 = cos(Phi);
    sp1 = sin(Phi);

    ct2 = cos(telescopeTheta);
    st2 = sin(telescopeTheta);
    cp2 = cos(telescopePhi);
    sp2 = sin(telescopePhi);

    x1 = st1*cp1; y1 = st1*sp1; z1 = ct1;
    x2 = st2*cp2; y2 = st2*sp2; z2 = ct2;

    *t1 = (telescopeTheta-Theta);
    *t2 = (telescopePhi-Phi);
    
    return (  acos(x1*x2 + y1*y2 + z1*z2) );
  }
  
  inline void print ( void ) {
    float *ptr = (float *)this;
    for(int i=0; i<SIZE_OF_MCEVENTHEADER_2; ++i,++ptr)
      cerr << i << ':' << *ptr << '\n';
  }

};
// @endcode

#endif  // not defined MCEventHeader_2_Class
