#ifndef __ROOT_R_E_TAB_ALL__
#define __ROOT_R_E_TAB_ALL__

#include <ctype.h>
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <string>
#include <fstream>
#include <iostream>
 
#include "TROOT.h"
#include "TSystem.h"




// [ATIMA | PRAL | PRAL-PROJ | VEDA | ZIEGLER | SRIM | PRAL-ZIEGLER | PRAL-SRIM]

// construct seperate objects for different sources of the tables
// (e.g. RE_ATIMA for ATIMA tables etc)

// Usage of Range-Energy tables created by ATIMA
// tables created using the script /u/lukasik/AIX/atima/job_si0:
//
//   atima batch jobs
// energy 0.01,2000,-200,log 0 0
// save satan
// mode fth
// par 1 18 19 24
// units MeV/u  mikrom  MeV/mg/cm2
// job
// H1    ----> Si,200 Si_H_1
// He4   ----> Si,200 Si_He_4
//        .........
// end job
//
// See az.dat for a full list of species
//
// CAUTION
// az.dat (or equivalent) is the "pilot" file containing species
// for which the tables were created in tab_dir
// There MUST be 1 to 1 correspondence between the species in
// the az.dat file and the tables stored in tab_dir 
//
// If the requested Ap and Zp is not present in the tables a search
// for the nearest species is performed (resulting in Atab and Ztab)
// and the following scaling of range is applied (see W.R. Leo p. 33):
/* 
                Ap  Ztab^2       /   Atab \
      Rp(Ep) = ---- ------ * Rtab|Ep ---- |
               Atab  Zp^2        \    Ap  /
*/ 
//
// NOTE
// for most abundant isotope use A=0;
// for Electron use Z=1, A=-1
// for Electron use Z=1, A=-2
// for Pion     use Z=1, A=-3
// for Muon     use Z=1, A=-4
// for Kaon     use Z=1, A=-5
//
// if(iz==1 && ia==-1) a=  0.511/931.49386;   //Gamma
// if(iz==1 && ia==-2) a=  0.511/931.49386;   //Electron
// if(iz==1 && ia==-3) a=139.56995/931.49386; //Pion
// if(iz==1 && ia==-4) a=105.65841/931.49386; //Muon
// if(iz==1 && ia==-5) a=493.67700/931.49386; //Kaon

//_____________________________________________________________________________ 

class TABLE{
private:
  static const int DIM=220;   // max number of points in the table
  float E[DIM], R[DIM], S[DIM];
  int   NSPMAXI;              // actual number of points in the table
  void  DERIVXY();
  void  READ_TABLE(       char *tab_fname, 
                    const char *tab_source,
                          char *spec, 
                           int  Z, 
                           int  A, 
                         float  M);
public:
  TABLE( char *tab_fname,
   const char *tab_source,
         char *spec, 
          int  Z, 
          int  A, 
        float  M);
  virtual ~TABLE();
  float  RANGE(float energy);
  float  ENERGY(float range);
  int    GetNSPMAXI(){return NSPMAXI;}
  float  GetR(int i){return (i<NSPMAXI && i>=0) ? R[i] : -1.;}
  float  GetE(int i){return (i<NSPMAXI && i>=0) ? E[i] : -1.;}
  float* GetR(){return &R[0];}
  float* GetE(){return &E[0];}

  ClassDef(TABLE,0) // Table 
};	        
//_____________________________________________________________________________ 

class RE_TAB{


private:
  static int ZA[92];                  // default mass number
  static const int NAZ=110;           // max number of species (tables) 
  int   NAZMAXI;                      // actual number of tables in use
  int   NAZMAXI_AZ;                   // number of species in az.dat
  int   A[NAZ],Z[NAZ];                // possible projectile A & Z (az.dat)
  float M[NAZ];                       // possible projectile mass (amu)
  char  SYM[NAZ][3];         
  const char *g_mat_prefix;           // material prefix (eg Si for Silicon)
  const char *g_tab_source;           // tables' source [ATIMA|PRAL...] 
  const char *g_tab_dir;              // tables' directory 
  TString     g_tab_dir_default_path; // path to calibration tables taken form environment varibale
  const char *g_az_spec;              // full name of the 'pilot' az.dat-like file
  void  INIT_APZP();
  int   NAPZP(int Zp,int Ap);         // returns consecutive number of apropriate
                                      // table[naz] for requested Zp, Ap
                            

  static const float AMU=931.4943;    // amu
  static const float MEl=0.511;       // Electron mass (MeV), use Z=1, A=-1
  static const float MPi=139.56995;   // Pion mass (MeV), use Z=1, A=-3
  static const float MMu=105.65841;   // Muon mass (MeV), use Z=1, A=-4
  static const float MK =493.677;     // Kaon mass (MeV), use Z=1, A=-5

  float MPART[5];                     // the above end up here

public:
  RE_TAB(const char *mat_prefix = "Si",
         const char *tab_source = "ATIMA",
         const char *tab_dir    = "default");
  virtual ~RE_TAB();

  struct TAB_STRUCT{
    TABLE *AZ;              // pointer to TABLE object
    int   Z;
    int   A;
    float M;
    char  SYM[3];
  };  
  TAB_STRUCT TAB[NAZ];
  
  float RANGE(int Zp, int Ap, float energy);
  float ENERGY(int Zp, int Ap, float range);
  void  INIT_TABLE(int naz);
  int	 GetNSPMAXI(int Zp, int Ap);
  const char*  GetSource(){return g_tab_source;}
  float  GetR(int Zp, int Ap, int i);
  float  GetE(int Zp, int Ap, int i);
  float* GetR(int Zp, int Ap);
  float* GetE(int Zp, int Ap);
  int    GetA(int Zp, int Ap);
  int    GetZ(int Zp, int Ap);
  float  GetM(int Zp, int Ap);
  char*  GetSym(int Zp, int Ap);
  int    A_EPAX(int Zp, int Zt, int At);
  int    A_GREEN(int Zp){
                        int a_gre;
                        a_gre = int((Zp-100.+sqrt(10000.+40.*Zp+Zp*Zp))/0.6);
			if(Zp==1) a_gre = 1;
			if(Zp==4) a_gre = 9;
                        return a_gre;}
  ClassDef(RE_TAB,0) // Range-Energy tables from ATIMA + scaling
};	        
//_____________________________________________________________________________ 

#endif