#ifndef PARTICLEBEAM_SEEN
#define  PARTICLEBEAM_SEEN
#include <iostream>
#include <sstream>
#include <stdio.h>
#include<cmath>
#include <vector>

#include "bareParticle.h"
#include "nomdElements.h"
#include "beam2Moments.h"

using namespace std;

class particleBeam 
{
  bool particleRepresentationOk_;
  bool momentRepresentationOk_;
  beam2Moments rij_;
  double P0Transport_;

  // ce faisceau contient les coordonnees des particules quand elles 
  // passe en un z, qui est celui de la particules de reference (
  // (particule donne par l'attribut  referenceParticle_ ci-dessous ; 
  // cette particule est abstraite, car elle est dupliquee dans le faisceau 
  // proprement dit, et n'intervient plus, en tant que telle, dans la dynamique)

  // Une particule courante est definie par : 
  // x (cm)
  // y (cm)
  // cdeltat = dephasage, donc décalage en z, exprime en c.dt soit, en cm  
  // cdeltat > 0 implique que la particule est en avance sur la particule de ref.
  // gamma.betax (radians ? a verifier)
  // gamma.betay
  // gamma.betaz

  vector<bareParticle> relativePartic_;
  // La particule de reference (abstraite) est definie par x,y,z,gamma.betax, gamma.betay, gamma.betaz
  // elle est dupliquee dans le faisceau proprement dit (relativePartic_) avec z = 0;

  bareParticle referenceParticle_;
  

  vector<double> centroid_;

  void impressionDesMoments() const;
  void razDesMoments();
  void particlesPhaseSpaceComponent(vector<double>& coord, unsigned index);
  unsigned indexFromPspaToTransport(unsigned index) const;

  void histogramInitialize(unsigned int index,vector<double>& vshf,double out[3]);
  void histogramPacking(double ecatyp,vector<double>vshf,vector<double>&xcor,vector<int>& hist);

  double getSigmaTransportij(unsigned indexI, unsigned indexJ);
double getUnnormalizedTranspPhaseSpaceArea(unsigned TranspIndexAbs, unsigned TranspIndexOrd);


 unsigned pspaCoorIndexFromString(string nameCoor) const;
 unsigned transportCoorIndexFromString(string nameCoor) const;

inline double dimensionalFactorFromTransportToGraphics(unsigned index) {
  if ( index == 0 || index == 2 || index == 4 ) return 10.;   // mm
  else return  1.0;  // mrad
}




 inline string transportVariableName(unsigned index) {
  switch ( index ) {
  case 0 : return  " x ";
  case 1 : return  " x' ";
  case 2 : return " y ";
  case 3 : return " y' ";
  case 4 : return " l ";
  case 5 : return "dp/p";
  default : { 
    cout << " particleBeam::transportVariableName : coordinate index ERROR index :  "<< index << endl;
    return "";
  }
  }
 }




 inline string graphicTransportUnitName(unsigned index) {
  switch ( index ) {
  case 0 : return  " mm ";
  case 1 : return  " mrad ";
  case 2 : return " mm ";
  case 3 : return " mrad ";
  case 4 : return " mm ";
  case 5 : return " % ";
  default : { 
    cout << " particleBeam::transportVariableName : coordinate index ERROR index :  "<< index << endl;
    return "";
  }
  }
 }


 public:

  particleBeam(); 
  virtual ~particleBeam() {;}

  void buildMomentRepresentation();
  void clear();
  int getNbParticles() const;
  const beam2Moments&  getTransportMoments() const;
  double getUnnormalizedEmittanceX();
  double getP0Transport() const;
  double referenceKineticEnergyMeV() const;
  bool particleRepresentationOk() const;
  bool momentRepresentationOk() const;
  void  addParticle( bareParticle p);
  const vector<bareParticle>& getParticleVector() const;
  vector<bareParticle>& getParticleVector(); 
  double getXmaxRms();
  //  void getVariance(double& varx, double& vary, double& varz) const;
  void set2Moments(beam2Moments& moments);
  void setWithParticles(vector<double>& centroid, bareParticle& referencePart, vector<bareParticle>& particles);
  void printAllXYZ() const;
  void ZrangeCdt(double& zmin, double& zmax) const; // extension en phase (cm)
  void donneesDessinEllipse(vector<double>& xcor, vector<double>& ycor, vector<string>& legende, string namex, string namey);
  void particlesPhaseSpaceData(vector<double>& xcor, vector<double>& ycor, vector<string>& legende, string namex, string namey);


 // coordonnees d'une particule dans le faisceau deploye ( passage 
 // d'une representation z=cte a une representation t=cte)
 TRIDVECTOR coordonneesDeployees(unsigned particleIndex, double cdtShift = 0.0);

 inline TRIDVECTOR betaGamma(unsigned particleIndex) {
     return relativePartic_.at(particleIndex).getBetaGamma();
   }

  virtual string fileOutputFlow() const;
  virtual bool FileInput(ifstream& ifs);
  virtual void writeToAMLFile(string fileName);
  virtual bool readFromAMLFile(string fileName);

  void histogramme(unsigned int index,vector<double>&xcor,vector<int>& hist,vector<string>& legende);
};
#endif