source: Sophya/trunk/SophyaPI/PIext/servnobjm.cc@ 333

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

#include <typeinfo>
#include <iostream.h>
#include <string>
#include <list>
#include <map>
#if defined(__KCC__)
using std::string ;
#include <list.h>
#include <map.h>
#endif

#include "strutil.h"

#include "nobjmgr.h"
#include "servnobjm.h"
#include "nomgadapter.h"
#include "pistdimgapp.h"

#include "fct1dfit.h"
#include "fct2dfit.h"

#include "matrix.h"
#include "cvector.h"
#include "ntuple.h"
#include "cimage.h"

#include "histos.h"
#include "histos2.h"
#include "ntuple.h"
#include "hisprof.h"

#include "piscdrawwdg.h"
#include "pisurfdr.h"

#include "pintuple.h"
#include "pintup3d.h"

#include "pipodrw.h"



/* --Methode-- */
Services2NObjMgr::Services2NObjMgr(NamedObjMgr* omg, string& tmpdir)
{
TmpDir = tmpdir;
PDynLinkMgr::SetTmpDir(tmpdir);
mImgapp = NULL;
mOmg = omg;
dynlink = NULL;
InitGrAttNames();
}
  
/* --Methode-- */
Services2NObjMgr::~Services2NObjMgr()
{
CloseDLL();
}

/* --Methode-- */
void Services2NObjMgr::RegisterClass(AnyDataObj* o, NObjMgrAdapter* oa)
{
ObjAdaptList::iterator it;
for(it = objadaplist.begin(); it != objadaplist.end(); it++) 
    if (typeid(*o) == typeid(*((*it).obj)))  THROW(dupIdErr);

dataobj_adapter oba;
oba.obj = o;
oba.obja = oa;
objadaplist.push_back(oba);
}

/* --Methode-- */
NObjMgrAdapter* Services2NObjMgr::GetAdapter(AnyDataObj* o)
{
ObjAdaptList::iterator it;
for(it = objadaplist.begin(); it != objadaplist.end(); it++) 
  if (typeid(*o) == typeid(*((*it).obj)))  return((*it).obja->Clone(o));
return(new NObjMgrAdapter(o));
}

/* --Methode-- */
void Services2NObjMgr::PlotFunc(string const & expfunc, string & nom, double xmin, double xmax, int np, string dopt)
{
FILE *fip;
string fname = TmpDir + "func1_pia_dl.c";
string  cmd;
int rc;

if (!mImgapp)  return; 

cmd = "rm -f " + fname;
rc = system(cmd.c_str());
// printf("PlotFunc_Do> %s  (Rc=%d)\n", cmd.c_str(), rc); 

if ((fip = fopen(fname.c_str(), "w")) == NULL)   { 
  string sn = fname; 
  cout << "Services2NObjMgr/PlotFunc_Erreur: Pb. Ouverture " << sn << endl; 
  return; 
  }

// constitution du fichier a compiler
fputs("#include <math.h> \n", fip);
fputs("double func1_pia_dl_func(double x) \n{\n", fip);
fprintf(fip,"return(%s); \n}\n", expfunc.c_str());
fclose(fip);

string func = "func1_pia_dl_func";
DlFunctionOfX f = (DlFunctionOfX) LinkFunctionFromFile(fname, func);
if (!f) return; 
PlotFunc(f, nom, xmin, xmax, np, dopt);
CloseDLL();
return;
}

/* --Methode-- */
void Services2NObjMgr::PlotFunc2D(string const & expfunc, string & nom, double xmin, double xmax, 
                                  double ymin, double ymax, int npx, int npy, string dopt)
{
FILE *fip;
string fname = TmpDir + "func2_pia_dl.c";
string cmd;
int rc;

if (!mImgapp)  return; 

cmd = "rm " + fname;
rc = system(cmd.c_str());
// printf("PlotFunc2D_Do> %s  (Rc=%d)\n", cmd.c_str(), rc); 

if ((fip = fopen(fname.c_str(), "w")) == NULL)   { 
  string sn = fname; 
  cout << "Services2NObjMgr/PlotFunc2D_Erreur: Pb. Ouverture " << sn << endl; 
  return; 
  }

// constitution du fichier a compiler
fputs("#include <math.h> \n", fip);
fputs("double func2_pia_dl_func(double x, double y) \n{\n", fip);
fprintf(fip,"return(%s); \n}\n", expfunc.c_str());
fclose(fip);

string func = "func2_pia_dl_func";
DlFunctionOfXY f = (DlFunctionOfXY) LinkFunctionFromFile(fname, func);
if (!f)  return;
PlotFunc2D(f, nom, xmin, xmax, ymin, ymax, npx, npy, dopt);
CloseDLL();
return;
}

/* --Methode-- */
void Services2NObjMgr::PlotFuncFrCFile(string const & fname, string const & func, string & nom, 
                                       double xmin, double xmax, int np, string dopt)
{
DlFunctionOfX f = (DlFunctionOfX) LinkFunctionFromFile(fname, func);
if (!f) return; 
PlotFunc(f, nom, xmin, xmax, np, dopt);
CloseDLL();
return;
}

/* --Methode-- */
void Services2NObjMgr::PlotFunc2DFrCFile(string const & fname, string const & func, string & nom, 
                                         double xmin, double xmax, double ymin, double ymax, int npx, int npy, string dopt)
{
DlFunctionOfXY f = (DlFunctionOfXY) LinkFunctionFromFile(fname, func);
if (!f)  return;
PlotFunc2D(f, nom, xmin, xmax, ymin, ymax, npx, npy, dopt);
CloseDLL();
return;
}

/* --Methode-- */
void Services2NObjMgr::PlotFunc(DlFunctionOfX f, string & nom, double xmin, double xmax, int np, string dopt)
{
if (!mImgapp)  return; 

int k;
if (np < 1) np = 1;
if (xmax <= xmin) xmax = xmin+1.;
Vector* vpy = new Vector(np);

double xx;
double dx = (xmax-xmin)/np;
TRY {
  for(k=0; k<np; k++) { xx = xmin+dx*k; (*vpy)(k) = f(xx); }
} CATCH(merr) {
  fflush(stdout); 
  cout << endl; 
  cerr << endl;
  string es = PeidaExc(merr);
  cerr << "Services2NObjMgr::PlotFunc()  Exception :" << merr << es;
  delete vpy;
  vpy = NULL;
  } ENDTRY;


if (vpy) {
  string titre = (nom.length() < 1) ? "Function f(x)" : nom;
  P1DArrayAdapter* vya = new POVectorAdapter(vpy, true);
  vya->DefineXCoordinate(xmin, (xmax-xmin)/np);
  PIYfXDrawer* dr = new   PIYfXDrawer(vya, NULL, true) ; 
  bool fgsr = true;
  dopt = "thinline," + dopt;
  int opt = DecodeDispOption(dopt, fgsr);
  int rsid = mImgapp->DispScDrawer(dr, titre, opt);
  if (fgsr) mImgapp->RestoreGraphicAtt(); 
  if (nom.length() > 0) { 
    mOmg->AddObj(vpy, nom);
    mOmg->AddWRsId(nom, rsid);
    }
  }

return;
}

/* --Methode-- */
void Services2NObjMgr::PlotFunc2D(DlFunctionOfXY f, string & nom, double xmin, double xmax, double ymin, double ymax, 
                             int npx, int npy, string dopt)
{
if (!mImgapp)  return; 

if (npx < 3) npx = 3;
if (npy < 3) npy = 3;
if (npx > 250) npx = 250;
if (npy > 250) npy = 250;
if (xmax <= xmin) xmax = xmin+1.;
if (ymax <= ymin) ymax = ymin+1.;

Matrix* mtx = new Matrix(npy, npx);

int i,j;
double xx, yy;
double dx = (xmax-xmin)/npx;
double dy = (ymax-ymin)/npy;
// printf(" -- DBG -- %d %d , %g %g , %g %g \n", npx, npy, xmin, xmax, ymin, ymax);
TRY {
  for(j=0; j<npy; j++) {
    yy = ymin+dy*j; 
    for(i=0; i<npx; i++) {
      xx = xmin+dx*i;
      (*mtx)(j, i) = f(xx, yy);
    }
  }
} CATCH(merr) {
  fflush(stdout); 
  cout << endl; 
  cerr << endl;
  string es = PeidaExc(merr);
  cerr << "Services2NObjMgr::PlotFunc2D()  Exception :" << merr << es;
  delete mtx;  mtx = NULL;
  } ENDTRY;

if (mtx) {
  string titre = (nom.length() < 1) ? "Function f(x,y)" : nom;  
  bool fgsr = true;
  int opt = DecodeDispOption(dopt, fgsr);
  P2DArrayAdapter* arr = new POMatrixAdapter(mtx, true);
  arr->DefineXYCoordinates(xmin, ymin, dx, dy);
  PISurfaceDrawer* sdr = new PISurfaceDrawer(arr, true, true, true);
  int rsid = mImgapp->Disp3DDrawer(sdr, titre, opt);
  if (fgsr) mImgapp->RestoreGraphicAtt();
  if (nom.length() > 0)  { 
    mOmg->AddObj(mtx, nom);
    mOmg->AddWRsId(nom, rsid);
    }
  }

return;
}

/* --Methode-- */
void Services2NObjMgr::DisplayPoints2D(string& nom,  string& expx, string& expy,
                                       string& experrx, string& experry,
                                       string& expcut, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::DisplayPoints2D() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

// Creation NTuple 
char* ntn[4] = {"expx","expy","expex","expey",};
NTuple* nt = NULL;
bool haserr = false;

if ( (experrx.length() > 0 ) && (experry.length() > 0 ) ) {  haserr = true;   nt = new NTuple(4, ntn); }
else {  haserr = false;  experrx = experry = "0."; nt = new NTuple(2, ntn); }

ComputeExpressions(obja, expx, expy, experrx, experry, expcut, nt, NULL, NULL);

if (nt->NEntry() < 1) {
  cout << "Services2NObjMgr::DisplayPoints2D() Warning  Zero points satisfy cut !" << endl; 
  delete nt;
  return;
  }

// nt->Show();
// nt->Print(0,10);
PINTuple* pin = new PINTuple(nt, true);
pin->SelectXY(ntn[0], ntn[1]);
if ( haserr ) pin->SelectErrBar(ntn[2], ntn[3]);

bool fgsr = true;
dopt = "defline," + dopt;
int opt = DecodeDispOption(dopt, fgsr);
string titre =  nom + ":" + expy + "%" + expx;
mImgapp->DispScDrawer( (PIDrawer*)pin, titre, opt);
if (fgsr) mImgapp->RestoreGraphicAtt();
return;
}

/* --Methode-- */
void Services2NObjMgr::DisplayPoints3D(string& nom,  string& expx, string& expy, 
                                       string& expz, string& expcut, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::DisplayPoints3D() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

char* ntn[3] = {"expx","expy","expz"};
NTuple* nt = new NTuple(3,ntn);  // Creation NTuple

string expwt = "1.";
ComputeExpressions(obja, expx, expy, expz, expwt, expcut, nt, NULL, NULL);

if (nt->NEntry() < 1) {
  cout << "Services2NObjMgr::DisplayPoints3D() Warning  Zero points satisfy cut !" << endl;
  delete nt;
  return;
  }
nt->Show();
nt->Print(0,10);
PINTuple3D* pin = new PINTuple3D(nt, true);
pin->SelectXYZ(ntn[0], ntn[1], ntn[2]);
bool fgsr = true;
dopt = "defline," + dopt;
int opt = DecodeDispOption(dopt, fgsr);

// Pour plot a partir de DispScDrawer
// string nomdisp = "_NT3D_";
// mImgapp->DispScDrawer( (PIDrawer*)pin, nomdisp, opt);
// Pour plot a partir de Disp3DDrawer
string titre =  nom + ":" + expy + "%" + expx;
mImgapp->Disp3DDrawer(pin, titre, opt);

if (fgsr) mImgapp->RestoreGraphicAtt();
return;
}

/* --Methode-- */
void Services2NObjMgr::DisplayPoints2DW(string& nom, string& expx, string& expy, 
                                        string& expwt, string& expcut, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::DisplayPoints2DW() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

char* ntn[3] = {"expx","expy","expw"};
NTuple* nt = new NTuple(3,ntn);  // Creation NTuple

string exp = "1.";
ComputeExpressions(obja, expx, expy, expwt, exp, expcut, nt, NULL, NULL);

if (nt->NEntry() < 1) {
  cout << "Services2NObjMgr::DisplayPoints2DW() Warning  Zero points satisfy cut !" << endl;
  delete nt;
  return;
  }

PINTuple* pin = new PINTuple(nt, true);
pin->SelectXY(ntn[0], ntn[1]);
pin->SelectWt(ntn[2]);

bool fgsr = true;
int opt = DecodeDispOption(dopt, fgsr);
string titre =  nom + ":" + expwt + "_" + expy + "%" + expx ;
mImgapp->DispScDrawer( (PIDrawer*)pin, titre, opt);
if (fgsr) mImgapp->RestoreGraphicAtt();
return;
}

/* --Methode-- */
void Services2NObjMgr::ProjectH1(string& nom, string& expx, string& expwt, 
                                 string& expcut, string& nomh1, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::ProjectH1() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

Histo* h1 = NULL;
NTuple* nt = NULL;
AnyDataObj* oh = NULL; 
if (nomh1.length() > 0) oh=mOmg->GetObj(nomh1);
else nomh1 = "/tmp/projh1d";
if ( (oh != NULL) && (typeid(*oh) == typeid(Histo)) )  h1 = (Histo*)oh;  // Pas de remise a zero ! h1->Zero();
else {
  char* ntn[2]= {"hxval", "hwt"};
  nt = new NTuple(2,ntn);  // Creation NTuple 
  }
string expz = "0.";   
ComputeExpressions(obja, expx, expwt, expz, expwt, expcut, nt, h1, NULL);

if ((!h1) && (!nt)) return;
if (!h1) {
  if (nt->NEntry() < 1) {
    cout << "Services2NObjMgr::ProjectH1() Warning  Zero points satisfy cut !" << endl; 
    delete nt;
    return;
    }
  double xmin, xmax;
  nt->GetMinMax(0, xmin, xmax);
  h1 = new Histo(xmin, xmax, 100);
  int k;
  float* xn;
  for(k=0; k<nt->NEntry(); k++)    {
    xn = nt->GetVec(k);
    h1->Add(xn[0], xn[1]);    
    }
  delete nt;
  mOmg->AddObj(h1, nomh1);
  }

mOmg->DisplayObj(nomh1, dopt);
return;  
}

/* --Methode-- */
void Services2NObjMgr::ProjectH2(string& nom, string& expx, string& expy, string& expwt, 
                                 string& expcut, string& nomh2, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::ProjectH2() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

Histo2D* h2 = NULL;
NTuple* nt = NULL;
AnyDataObj* oh = NULL;
if (nomh2.length() > 0)  oh=mOmg->GetObj(nomh2);
else nomh2 = "/tmp/projh2d";
if ( (oh != NULL) && (typeid(*oh) == typeid(Histo2D)) )  h2 = (Histo2D*)oh;  // Pas de remise a zero ! h2->Zero();
else {
  char* ntn[3]= {"hxval", "hyval", "hwt"};
  nt = new NTuple(3,ntn);  // Creation NTuple 
  }
string expz = "0.";   
ComputeExpressions(obja, expx, expy, expwt, expwt, expcut, nt, NULL, h2);

if ((!h2) && (!nt)) return;
if (!h2) {
  if (nt->NEntry() < 1) {
    cout << "Services2NObjMgr::ProjectH2() Warning  Zero points satisfy cut !" << endl; 
    delete nt;
    return;
    }
  double xmin, xmax, ymin, ymax;
  nt->GetMinMax(0, xmin, xmax);
  nt->GetMinMax(0, ymin, ymax);
  h2 = new Histo2D(xmin, xmax, 50, ymin, ymax, 50);
  int k;
  float* xn;
  for(k=0; k<nt->NEntry(); k++)    {
    xn = nt->GetVec(k);
    h2->Add(xn[0], xn[1], xn[2]);    
    }
  delete nt;
  mOmg->AddObj(h2, nomh2);
  }

mOmg->DisplayObj(nomh2, dopt);
return;  

}

/* --Methode-- cmv 13/10/98 */
void Services2NObjMgr::ProjectHProf(string& nom, string& expx, string& expy, string& expwt,
                                    string& expcut, string& nomprof, string dopt)
//      Pour remplir un ``GeneralFitData'' a partir de divers objets:
//| nom = nom de l'objet a projeter dans un HProf.
//| expx = expression X de definition du bin.
//| expy = expression Y a additionner dans le bin.
//| expwt = expression W du poids a additionner.
//| expcut = expression du test de selection.
//| nomprof = nom du HProf engendre (optionnel). Si l'objet n'existe pas
//|           les limites Xmin,Xmax sont calculees automatiquement.
//|           sinon ce sont celles de l'objet preexistant.
//| opt = options generales pour le display.
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::ProjectHProf() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

HProf* hprof = NULL;
NTuple* nt = NULL;
AnyDataObj* oh = NULL;
if (nomprof.length() > 0)  oh=mOmg->GetObj(nomprof);
else nomprof = "/tmp/projprof";
if( (oh!=NULL) && (typeid(*oh) == typeid(HProf)) )  hprof = (HProf*)oh;
else {
  char* ntn[3]= {"hxval", "hyval", "hwt"};
  nt = new NTuple(3,ntn);  // Creation NTuple 
}
string expz = "0.";
ComputeExpressions(obja, expx, expy, expwt, expwt, expcut, nt, NULL, NULL, hprof);

if((!hprof) && (!nt)) return;
if(!hprof) {
  if (nt->NEntry() < 1) {
    cout << "Services2NObjMgr::ProjectHProf() Warning  Zero points satisfy cut !" << endl;
    delete nt;
    return;
  }
  double xmin, xmax;
  nt->GetMinMax(0, xmin, xmax);
  hprof = new HProf(xmin, xmax, 100);
  int k;
  float* xn;
  for(k=0; k<nt->NEntry(); k++)    {
    xn = nt->GetVec(k);
    hprof->Add(xn[0], xn[1], xn[2]);
    }
  delete nt;
  mOmg->AddObj(hprof, nomprof);
  }
hprof->UpdateHisto();

mOmg->DisplayObj(nomprof, dopt);
return;
}

/* --Methode-- */
void Services2NObjMgr::FillVect(string& nom, string& expx, string& expcut, 
                                string& nomvec, string dopt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::FillVect() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

NTuple* nt = NULL;
if (nomvec.length() < 1) nomvec = "/tmp/fillvec";

char* ntn[2]= {"vecval", "vecwt"};
nt = new NTuple(1,ntn);  // Creation NTuple
 
string expwt = "1.";
string expz = "0.";   
ComputeExpressions(obja, expx, expz, expz, expwt, expcut, nt, NULL, NULL);

if (!nt) return;
if (nt->NEntry() < 1) {
  cout << "Services2NObjMgr::FillVect() Warning  Zero points satisfy cut !" << endl; 
  delete nt;
  return;
  }
  
Vector* vec = new Vector(nt->NEntry());
int k;
float* xn;
for(k=0; k<nt->NEntry(); k++)    {
  xn = nt->GetVec(k);
  (*vec)(k) = xn[0];
  }
delete nt;
mOmg->AddObj(vec, nomvec); 
mOmg->DisplayObj(nomvec, dopt);
return;  
}

/* --Methode-- */
void Services2NObjMgr::FillNT(string& nom, string& expx, string& expy, string& expz, 
                              string& expt, string& expcut, string& nomnt)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::FillNT() Error , Pas d'objet de nom " << nom << endl; 
  return;
  }
if (!mImgapp)  return; 

bool fgnnt = false;
NTuple* nt = NULL;
AnyDataObj* oh = NULL;
if (nomnt.length() > 0)  oh=mOmg->GetObj(nomnt);
else nomnt = "/tmp/fillnt";
if ( (oh != NULL) && (typeid(*oh) == typeid(NTuple)) )  { 
  nt = (NTuple*)oh;
  if (nt->NVar() > 10) {
    cout << "Services2NObjMgr::FillNT() Warning , Max 10 var ds NTuple -> new NTuple" << endl;
    nt = NULL;
    }
  }  
if (nt == NULL) {
  char* ntn[4]= {"x", "y","z","t"};
  nt = new NTuple(4,ntn);  // Creation NTuple 
  fgnnt = true;
  }

ComputeExpressions(obja, expx, expy, expz, expt, expcut, nt, NULL, NULL);

if (fgnnt) mOmg->AddObj(nt, nomnt);
return;  

}

/* --Methode-- */
void Services2NObjMgr::FillNTFrCFile(string & nom, string const & fname, 
                                     string const & funcname, string & nomnt)
{
if (!mImgapp)  return; 

NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::FillNTFrCFile( " << nom << "...) No such object" <<endl;
  return;
  }
NTupleInterface* objnt = obja->GetNTupleInterface();
if (objnt == NULL)  {
  cout << "Services2NObjMgr::FillNTFrCFile( " << nom << "...) No NTupleInterface !" <<endl;
  return;
  }

NTLoopExprFunc f = (NTLoopExprFunc)LinkFunctionFromFile(fname, funcname);
if (!f) { 
  cerr << "Services2NObjMgr::FillNTFrCFile Error Creation  NTLoopExprFunc" <<  endl; 
  return;
  }

bool fgnnt = false;
NTuple* nt = NULL;
if (nomnt.length() > 0) {
  AnyDataObj* oh = NULL;
  oh=mOmg->GetObj(nomnt);
  if ( (oh != NULL) && (typeid(*oh) == typeid(NTuple)) )  { 
    nt = (NTuple*)oh;
    if (nt->NVar() > 10) {
      cout << "Services2NObjMgr::FillNTFrCFile() Warning , Max 10 var ds NTuple -> new NTuple" << endl;
      nt = NULL;
      }
    }  
  if (nt == NULL) {
    char* ntn[4]= {"x", "y","z","t"};
    nt = new NTuple(4,ntn);  // Creation NTuple 
    fgnnt = true;
    }
  }

double xnt[10];
float  fxnt[10];

int i,k;
for(i=0; i<10; i++) fxnt[i] = xnt[i] = 0.;


// $CHECK$ A virer des que possible  - Pb blocage application quand trop d'impression
// redirige -    On redirige la sortie sur le terminal 
bool red = mImgapp->HasRedirectedStdOutErr();
mImgapp->RedirectStdOutErr(false);

TRY  {
  double* xn;
  int kmax = objnt->NbLines();
  for(k=0; k<kmax; k++)    {
    xn = objnt->GetLineD(k);
    if (f(xn, xnt, xnt+1, xnt+2, xnt+3, k, kmax) != 0) {
      if (nt) { 
         for(i=0; i<4; i++) fxnt[i] = xnt[i]; 
         nt->Fill(fxnt);
         }
      }
    }
 }
CATCH(merr) {
  fflush(stdout); 
  cout << endl; 
  cerr << endl;
  string es = PeidaExc(merr);
  cerr << "Services2NObjMgr::FillNTFrCFile()  Exception :" << merr << es;
  } ENDTRY;
 
CloseDLL();

// $CHECK$ A virer des que possible  On redirige la sortie sur la fenetre PIConsole 
mImgapp->RedirectStdOutErr(red);

if (fgnnt) mOmg->AddObj(nt, nomnt);
return;
}

/* --Methode-- */
void Services2NObjMgr::PrepareNTExpressionCFile(string & nom, string const & fname, 
                                                string const & funcname)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::PrepareNTExpressionCFile( " << nom << "...) No such object" <<endl;
  return;
  }
NTupleInterface* objnt = obja->GetNTupleInterface();
if (objnt == NULL)  {
  cout << "Services2NObjMgr::PrepareNTExpressionCFile( " << nom 
       << "...) No NTupleInterface !" <<endl;
  return;
  }
string  vardec = objnt->VarList_C("_xnti_");

FILE *fip;
if ((fip = fopen(fname.c_str(), "w")) == NULL)   { 
  cout << "Services2NObjMgr::PrepareNTExpressionCFile()_Error: fopen " << fname << endl; 
  return; 
  }
// constitution du fichier des decalarations des variables de l'interface NTuple
fputs("#include <stdlib.h> \n", fip);
fputs("#include <stdio.h> \n", fip);
fputs("#include <math.h> \n\n", fip);
fputs("/* NTupleInterface Variable declaration - Generated by piapp  \n", fip);
fputs("      -- Services2NObjMgr::PrepareNTExpressionCFile()  --     */ \n\n", fip);
fprintf(fip,"int %s(double* _xnti_, double* _rx_, double* _ry_, double* _rz_, \n", 
        funcname.c_str());
fprintf(fip,"       double* _rt_, int _n_, int _nmax_) \n");
fprintf(fip, "{ \n %s \n", vardec.c_str());
fputs("  if (!1) { /* Cut Expression failed */ \n", fip);
fputs("    *_rx_ = *_ry_ = *_rz_ = *_rt_ = 0.;  return(0);", fip);
fputs("  } \n  /* Cut expression satisfied */ \n", fip);
fputs("  *_rx_ = 1.; \n  *_ry_ = 1.; \n  *_rz_ = 1.; \n  *_rt_ = 1.; \n", fip);
fputs("  return(1); \n} \n", fip);

fclose(fip);
return;
}


/* --Methode-- cmv 13/10/98 */
void Services2NObjMgr::FillGFD(string& nom, string& expx, string& expy, string& expz,
                               string& experr, string& expcut, string& nomgfd)
//      Pour remplir un ``GeneralFitData'' a partir de divers objets:
//| nom = nom de l'objet a transcrire selon 1D: Z=f(X) ou 2D: Z=f(X,Y) .
//|       Vector,Matrix,Histo,HProf,Histo2D,Image<T>,StarList,NTuple,GeneralFitData
//| expx = expression X du GeneralFitData (1er abscisse)
//| expy = expression Y du GeneralFitData (2sd abscisse si non "", Z=f(X,Y))
//| expz = expression Z du GeneralFitData (valeur de l'ordonnee)
//| experr = expression de l'erreur sur l'ordonnee Z
//| expcut = expression du test de selection
//| nomgfd = nom du GeneralFitData engendre (optionnel)
{
NObjMgrAdapter* obja=NULL;
obja = mOmg->GetObjAdapter(nom);
if (obja == NULL) {
  cout << "Services2NObjMgr::FillGFD() Error , Pas d'objet de nom "<<nom<<endl;
  return;
  }
if(!mImgapp)  return;

// 2D ou 3D?
int nvar = 2;
if(expy.length()<=0) {nvar = 1; expy = "0.";}

// $CHECK$ - cmv calculait le nombre d'entree ndata
// en fonction de chaque objet - Je l'ai vire Reza 11/05/99

// Creation NTuple Buffer
char* ntn[4]= {"x","y","f","e"};
NTuple*nt = new NTuple(4,ntn);

// Remplissage NTuple buffer
ComputeExpressions(obja, expx, expy, expz, experr, expcut, nt, NULL, NULL);
if(nt->NEntry() < 1)
  {cout<<"Services2NObjMgr::FillGFD() Warning  Zero points satisfy cut !"<<endl; 
   delete nt; return;}

//Remplissage de la structure GeneraFitData
if (nt->NEntry() <= 0) {
  cout<<"Services2NObjMgr::FillGFD() Warning - NData= " << nt->NEntry() << endl;
  delete nt;
  return;
  }

GeneralFitData* gfd = new GeneralFitData(nvar,nt->NEntry(),0);
int k;
float* xn;
for(k=0; k<nt->NEntry(); k++) {
  xn = nt->GetVec(k);
  gfd->AddData(xn,xn[2],xn[3]);
}

// Menage et table d'objets
delete nt;
mOmg->AddObj(gfd, nomgfd);
return;
}


///////////////////// Fit 1D et 2D //////////////////////////
/* --Function static propre aux routines de fit 1D et 2D-- cmv 13/10/98 */
struct DFOptions {
  int okres, okfun;
  int polcx,polcy; double xc,yc;
  double err_e, err_E;
  double stc2;
  int lp,lpg;
  int i1,i2,j1,j2;
};
typedef struct DFOptions DFOPTIONS;
static void DecodeFitsOptions(string par,string step,string min,string max,string opt
                      ,Vector& Par,Vector& Step,Vector& Min,Vector& Max,DFOPTIONS& O);
void DecodeFitsOptions(string par,string step,string min,string max,string opt
                      ,Vector& Par,Vector& Step,Vector& Min,Vector& Max,DFOPTIONS& O)
//| Pour decoder les "string" et remplir les vecteurs du fit (cf commentaires dans Fit1D)
{
// set des vecteurs et decodage des string correspondantes
int NParMax = 100;
Par.Realloc(NParMax); Step.Realloc(NParMax);
Min.Realloc(NParMax); Max.Realloc(NParMax);
{
  Vector* v=NULL; string* s=NULL;
  {for(int i=0;i<NParMax;i++) {Par(i)=0.; Step(i)=1.; Min(i)=1.; Max(i)=-1.;}}
  for(int j=0;j<4;j++) {
    if(j==0)      {v=&Par; s=&par;}
    else if(j==1) {v=&Step; s=&step;}
    else if(j==2) {v=&Min; s=&min;}
    else if(j==3) {v=&Max; s=&max;}
    if(s->length()>0) *s += ",";
    for(int i=0;i<NParMax;i++) {
      if(s->length()<=0) break;
      sscanf(s->c_str(),"%lf",&(*v)(i));
      size_t p = s->find_first_of(',') + 1;
      if(p>=s->length()) *s = ""; else *s = s->substr(p);
    }
  }
}

// Decodage de options de opt
O.okres = O.okfun = 0;
O.polcx = O.polcy = 0;
O.xc = O.yc = 0.;
O.stc2 = 1.e-3;
O.err_e = O.err_E = -1.;
O.lp = 1; O.lpg = 0;
O.i1 = O.j1 = O.i2 = O.j2 = -1;

if(opt.length()<=0) return;
opt = "," + opt + ",";

if(strstr(opt.c_str(),",r,")) O.okres = 1;  // residus
if(strstr(opt.c_str(),",f,")) O.okfun = 1;  // fonction fittee
if(strstr(opt.c_str(),",x")) { // Demande de centrage (fit X=x-xc)
  O.polcx = 2; // Le centrage est calcule automatiquement
  size_t p = opt.find(",x");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) {
    sscanf(dum.c_str(),",x%lf",&O.xc);
    O.polcx = 1; // Le centrage est fixe par la valeur lue
  }
}
if(strstr(opt.c_str(),",y")) { // Demande de centrage (fit Y=y-yc)
  O.polcy = 2; // Le centrage est calcule automatiquement
  size_t p = opt.find(",y");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) {
    sscanf(dum.c_str(),",y%lf",&O.yc);
    O.polcy = 1; // Le centrage est fixe par la valeur lue
  }
}
if(strstr(opt.c_str(),",E")) { // Erreurs imposees a "sqrt(val)" ou "aa.b*sqrt(val)"
  size_t p = opt.find(",E");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) sscanf(dum.c_str(),",E%lf",&O.err_E);
  if(O.err_E<=0.) O.err_E = 1.;
  O.err_e=-1.;
}
if(strstr(opt.c_str(),",e")) { // Erreurs imposees a "1" ou "aa.b"
  size_t p = opt.find(",e");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) sscanf(dum.c_str(),",e%lf",&O.err_e);
  if(O.err_e<=0.) O.err_e = 1.;
  O.err_E=-1.;
}
if(strstr(opt.c_str(),",X")) { // Valeur du StopChi2
  size_t p = opt.find(",X");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) sscanf(dum.c_str(),",X%lf",&O.stc2);
  if(O.stc2<=0.) O.stc2 = 1.e-3;
}
if(strstr(opt.c_str(),",l")) { // niveau de print
  size_t p = opt.find(",l");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  float ab;
  if(dum.length()>2) sscanf(dum.c_str(),",l%f",&ab);
  if(ab<0) ab = 0.;
  O.lp = (int) ab; O.lpg = (int) 10.*(ab-O.lp);
}
if(strstr(opt.c_str(),",I")) { // intervalle de fit selon X
  size_t p = opt.find(",I");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) sscanf(dum.c_str(),",I%d/%d",&O.i1,&O.i2);
}
if(strstr(opt.c_str(),",J")) { // intervalle de fit selon Y
  size_t p = opt.find(",J");
  size_t q = opt.find_first_of(',',p+1);
  string dum = opt.substr(p,q-p);
  if(dum.length()>2) sscanf(dum.c_str(),",J%d/%d",&O.j1,&O.j2);
}
return;
}

/* --Methode-- cmv 13/10/98 */
void  Services2NObjMgr::Fit12D(string& nom, string& func,
                               string par,string step,string min,string max,
                               string opt)
//| --------------- Fit d'objets a 1 et 2 dimensions ---------------
//| nom  : nom de l'objet qui peut etre:
//|        fit-1D:  Vector,Histo1D,HProf ou GeneraFitData(1D)
//|        fit-2D:  Matrix,Histo2D,Imagexx ou GeneraFitData(2D)
//| func : pnn = fit polynome degre nn avec classe Poly (lineaire) 1D ou 2D
//|      : Pnn = fit polynome degre nn avec GeneralFit (non-lineaire) 1D ou 2D
//|      : gnn = fit gaussienne (hauteur) + polynome de degre nn 1D
//|      : g   = fit gaussienne (hauteur) 1D
//|      : enn = fit exponentielle + polynome de degre nn 1D
//|      : e   = fit exponentielle 1D
//|      : Gnn = fit gaussienne (volume) + polynome de degre nn 1D
//|      : G   = fit gaussienne (volume) 1D
//|      :     = fit gaussienne+fond (volume) 2D
//|      : Gi  = fit gaussienne+fond integree (volume) 2D
//|      : d   = fit DL de gaussienne+fond (volume) 2D
//|      : di  = fit DL de gaussienne+fond integree (volume) 2D
//|      : D   = fit DL de gaussienne+fond avec coeff variable p6 (volume) 2D
//|      : Di  = fit DL de gaussienne+fond integree avec coeff variable p6 (volume) 2D
//|      : M   = fit Moffat+fond (expos=p6) (volume) 2D
//|      : Mi  = fit Moffat+fond integree (expos=p6) (volume) 2D
//| par  : p1,...,pn valeur d'initialisation des parametres (def=0)
//| step : s1,...,sn valeur des steps de depart (def=1)
//| min  : m1,...,mn valeur des minima (def=1)
//| max  : M1,...,Mn valeur des maxima (def=-1) (max<=min : pas de limite)
//| opt  : options "Eaa.b,eaa.b,f,r,caa.b,Xaa.b"
//|       f = generation d'un Objet identique contenant la fonction fittee
//|       r = generation d'un Objet identique contenant les residus
//|       Xaa.b = aa.b valeur du DXi2 d'arret (def=1.e-3)
//|       la.b = niveau "a.b" de print: a=niveau de print Fit1/2D
//|                                     b=niveau de debug GeneralFit
//|       Ii1/i2 numeros des bins X de l'histos utilises pour le fit [i1,i2]
//|2D     Jj1/j2 numeros des bins Y de l'histos utilises pour le fit [j1,j2]
//|       - L'erreur est celle associee a l'objet si existe, 1 sinon sauf si
//|       E = erreur prise comme la racine de la valeur a fitter
//|       Eaa.b = erreur prise aa.b*sqrt(val)
//|       e = erreur prise egale a 1 pour toutes les valeurs
//|       eaa.b = erreur prise egale a aa.b
//|       xaa.b = demande de centrage: on fit x-aa.b au lieu de x)
//|       x = demande de centrage: on fit x-xc au lieu de x
//|           avec xc=abscisse du milieu de l'histogramme
//|           Actif pour: exp+poly 1D, poly 1D
//|                       gauss+poly 1D (mais xc est le centre de la gaussienne)
//|2D     yaa.b et y = idem "xaa.b et x" mais pour y
{
AnyDataObj* obj=mOmg->GetObj(nom);
if (obj == NULL) {
  cout<<"Services2NObjMgr::Fit12D() Error , Pas d'objet de nom "<<nom<<endl;
  return;
}
if (!mImgapp)  return;
if(func.length()<=0)
  {cout<<"Services2NObjMgr::Fit12D() Donnez un nom de fonction a fitter."<<endl;
   return;}
string ctyp = typeid(*obj).name();

int ndim = 0, nbinx=0, nbiny=0, ndata = 0;
Vector* v = NULL; Histo* h = NULL;
Matrix* m = NULL; Histo2D* h2 = NULL; RzImage* im = NULL;
GeneralFitData* g = NULL;

  // 1D
if (typeid(*obj) == typeid(Vector)) {
  ndim = 1;
  v = (Vector*) obj; nbinx = v->NElts(); nbiny = 1;
  }
else if ( (typeid(*obj) == typeid(HProf)) || (typeid(*obj) == typeid(Histo)) ) {
  ndim = 1;
  h = (Histo*)  obj; nbinx = h->NBins(); nbiny = 1;
  }
else if (typeid(*obj) == typeid(Matrix)) {
  ndim = 2;
  m = (Matrix*) obj; nbinx = m->NCol(); nbiny = m->NRows();
  }
else if (typeid(*obj) == typeid(Histo2D)) {
  ndim = 2;
  h2 = (Histo2D*) obj; nbinx = h2->NBinX(); nbiny = h2->NBinY();
  }
else if (typeid(*obj) == typeid(GeneralFitData)) {
  g = (GeneralFitData*) obj; nbinx = g->NData(); nbiny = 1;
  if(     g->NVar()==1) ndim = 1;
  else if(g->NVar()==2) ndim = 2;
  else {
     cout<<"GeneralFitData ne peut avoir que 1 ou 2 variables d'abscisse: "
         <<((GeneralFitData*) obj)->NVar()<<endl; return; }
  }
else if (dynamic_cast<RzImage*>(obj)) {
  ndim = 2;
  im = (RzImage*) obj; nbinx = im->XSize(); nbiny = im->YSize();
  }
else  {
  cout<<"Services2NObjMgr::Fit12D() Error , Objet n'est pas un "
      <<"Histo1D/HProf/Vector/Histo2D/Image/Matrix/GeneralFitData "<<ctyp<<endl;
  return; 
  }

ndata = nbinx*nbiny;
if(ndata<=0)
  {cout<<"L'objet a "<<nbinx<<","<<nbiny<<" bins ("<<ndata<<")"<<endl; return;}

// Decodage des options et des parametres, mise en forme
Vector Par(1); Vector Step(1); Vector Min(1); Vector Max(1); DFOPTIONS O;
DecodeFitsOptions(par,step,min,max,opt,Par,Step,Min,Max,O);
O.i1 = (O.i1<0||O.i1>=nbinx)? 0: O.i1;
O.i2 = (O.i2<0||O.i2>=nbinx||O.i2<O.i1)? nbinx-1: O.i2;
if(ndim>=2) {
  O.j1 = (O.j1<0||O.j1>=nbiny)? 0: O.j1;
  O.j2 = (O.j2<0||O.j2>=nbiny||O.j2<O.j1)? nbiny-1: O.j2;
} else O.j2 = O.j1 = 0;
if(O.polcx==2) {
  if(v||m)    O.xc = (O.i2-O.i1+1)/2.;
  else if(h)  O.xc = (h->XMin()+h->XMax())/2.;
  else if(h2) O.xc = (h2->XMin()+h2->XMax())/2.;
  else if(g)  {double mini,maxi; g->GetMinMax(2,mini,maxi); O.xc=(mini+maxi)/2.;}
  else if(im) {O.xc = im->XOrg() * im->XPxSize()*(O.i2-O.i1+1)/2.;}
}
if(O.polcy==2 && ndim>=2) {
  if(m)  O.yc = (O.j2-O.j1+1)/2.;
  if(h2) O.yc = (h2->YMin()+h2->YMax())/2.;
  if(g)  {double mini,maxi; g->GetMinMax(12,mini,maxi); O.yc=(mini+maxi)/2.;}
  if(im) {O.yc = im->YOrg() * im->YPxSize()*(O.j2-O.j1+1)/2.;}
}
if(O.lp>0)
  cout<<"Fit["<<nbinx<<","<<nbiny<<"] ("<<ndata<<") dim="<<ndim<<":"
      <<" Int=["<<O.i1<<","<<O.i2<<"],["<<O.j1<<","<<O.j2<<"]"<<endl
      <<" Cent="<<O.polcx<<","<<O.polcy<<","<<O.xc<<"+x"<<","<<O.yc<<"+y"
      <<" TypE="<<O.err_e<<","<<O.err_E
      <<" StpX2="<<O.stc2
      <<" lp,lpg="<<O.lp<<","<<O.lpg<<endl;

///////////////////////////////////
// Remplissage de GeneralFitData //
///////////////////////////////////
GeneralFitData mydata(ndim,ndata,0);
{for(int i=O.i1;i<=O.i2;i++) for(int j=O.j1;j<=O.j2;j++) {
  double x,y,f,e;

  if(v)
    {x= (double) i; f=(*v)(i); e=1.;}
  else if(h)
    {x=h->BinCenter(i); f=(*h)(i); e=(h->HasErrors())?h->Error(i):1.;}
  else if(m)
    {x=(double) i; y=(double) j; f=(*m)(j,i); e=1.;}
  else if(h2)
    {float xf,yf; h2->BinCenter(i,j,xf,yf); x=(double)xf; y=(double)yf;
     f=(*h2)(i,j); e=(h2->HasErrors())?h2->Error(i,j):1.;}
  else if(im)
    {x=im->XOrg()+(i+0.5)*im->XPxSize(); y=im->YOrg()+(j+0.5)*im->YPxSize();
     f=im->DValue(i,j); e=1.;}
  else if(g&&ndim==1) {x= g->X(i); f=g->Val(i); e=g->EVal(i);}
  else if(g&&ndim==2) {x= g->X(i); y= g->Y(i); f=g->Val(i); e=g->EVal(i);}
  else x=y=f=e=0.;

  // Gestion des erreurs a utiliser
  if(O.err_e>0.) e=O.err_e;
  else if(O.err_E>0.) {e=(y<-1.||y>1.)?O.err_E*sqrt(fabs(y)):O.err_E;}

  // Remplissage de generalfit
  if(func[0]=='p') {x -= O.xc; if(ndim>=2) y -= O.yc;}
  if(ndim==1)      mydata.AddData1(x,f,e);
  else if(ndim==2) mydata.AddData2(x,y,f,e);
}}
if(mydata.NData()<=0)
  {cout<<"Pas de donnees dans GeneralFitData: "<<mydata.NData()<<endl;
   return;}
if(O.lpg>1) {
  mydata.PrintStatus();
  mydata.PrintData(0);
  mydata.PrintData(mydata.NData()-1);
}

////////////////////////////////////////////
// Identification de la fonction a fitter //
////////////////////////////////////////////
GeneralFunction* myfunc = NULL;
if(func[0]=='p' && ndim==1) {
  // Fit de polynome sans passer par les GeneralFit
  int degre = 0;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit (lineaire) 1D polynome de degre "<<degre<<endl;
  Poly p1(0);
  double c2rl = mydata.PolFit(0,p1,degre);
  cout<<"C2r_lineaire = "<<c2rl<<endl;
  if(O.lp>0) cout<<p1<<endl;
  return;

} else if(func[0]=='P' && ndim==1) {
  // Fit de polynome
  int degre = 0;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit polynome 1D de degre "<<degre<<endl;
  Polyn1D* myf = new Polyn1D(degre,O.xc);
  myfunc = myf;

} else if(func[0]=='e' && ndim==1) {
  // Fit d'exponentielle
  int degre =-1;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit d'exponentielle+polynome 1D de degre "<<degre<<endl;
  Exp1DPol* myf;
  if(degre>=0) myf = new Exp1DPol((unsigned int)degre,O.xc);
       else    myf = new Exp1DPol(O.xc);
  myfunc = myf;

} else if(func[0]=='g' && ndim==1) {
  // Fit de gaussienne en hauteur
  int degre =-1;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit de Gaussienne_en_hauteur+polynome 1D de degre "<<degre<<endl;
  Gauss1DPol* myf;
  if(degre>=0) myf = new Gauss1DPol((unsigned int)degre,((O.polcx)?true:false));
  else { bool bfg = (O.polcx)?true:false;   myf = new Gauss1DPol(bfg); }
  myfunc = myf;

} else if(func[0]=='G' && ndim==1) {
  // Fit de gaussienne en volume
  int degre =-1;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit de Gaussienne_en_volume+polynome 1D de degre "<<degre<<endl;
  GaussN1DPol* myf;
  if(degre>=0) myf = new GaussN1DPol((unsigned int)degre,((O.polcx)?true:false));
  else  { bool bfg = (O.polcx)?true:false;   myf = new GaussN1DPol(bfg); }
  myfunc = myf;

} else if(func[0]=='p' && ndim==2) {
  // Fit de polynome 2D sans passer par les GeneralFit
  int degre = 0;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit (lineaire) polynome 2D de degre "<<degre<<endl;
  Poly2 p2(0);
  double c2rl = mydata.PolFit(0,1,p2,degre);
  cout<<"C2r_lineaire = "<<c2rl<<endl;
  if(O.lp>0) cout<<p2<<endl;
  return;

} else if(func[0]=='P' && ndim==2) {
  // Fit de polynome 2D
  int degre = 0;
  if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
  cout<<"Fit polynome 2D de degre "<<degre<<endl;
  Polyn2D* myf = new Polyn2D(degre,O.xc,O.yc);
  myfunc = myf;

} else if(func[0]=='G' && ndim==2) {
  // Fit de gaussienne+fond en volume
  int integ = 0;
  if(func.length()>1) if(func[1]=='i') integ=1;
  cout<<"Fit de Gaussienne+Fond 2D integ="<<integ<<endl;
  if(integ) {GauRhInt2D* myf = new GauRhInt2D; myfunc = myf;}
    else    {GauRho2D* myf = new GauRho2D; myfunc = myf;}

} else if(func[0]=='d' && ndim==2) {
  // Fit de DL gaussienne+fond en volume
  int integ = 0;
  if(func.length()>1) if(func[1]=='i') integ=1;
  cout<<"Fit de DL de Gaussienne+Fond 2D integ="<<integ<<endl;
  if(integ) {GdlRhInt2D* myf = new GdlRhInt2D; myfunc = myf;}
    else    {GdlRho2D* myf = new GdlRho2D; myfunc = myf;}

} else if(func[0]=='D' && ndim==2) {
  // Fit de DL gaussienne+fond avec coeff variable p6 en volume
  int integ = 0;
  if(func.length()>1) if(func[1]=='i') integ=1;
  cout<<"Fit de DL de Gaussienne+Fond avec coeff variable (p6) 2D integ="<<integ<<endl;
  if(integ) {Gdl1RhInt2D* myf = new Gdl1RhInt2D; myfunc = myf;}
    else    {Gdl1Rho2D* myf = new Gdl1Rho2D; myfunc = myf;}

} else if(func[0]=='M' && ndim==2) {
  // Fit de Moffat+fond (volume)
  int integ = 0;
  if(func.length()>1) if(func[1]=='i') integ=1;
  cout<<"Fit de Moffat+Fond (expos=p6) 2D integ="<<integ<<endl;
  if(integ) {MofRhInt2D* myf = new MofRhInt2D; myfunc = myf;}
    else    {MofRho2D* myf = new MofRho2D; myfunc = myf;}

} else {
  cout<<"Fonction "<<func<<" inconnue pour la dim "<<ndim<<endl;
  return;
}

/////////////////////////
// Fit avec generalfit //
/////////////////////////
if(myfunc->NPar()>Par.NElts())
  {cout<<"Trop de parametres: "<<myfunc->NPar()<<">"<<Par.NElts()<<endl;
  if(myfunc) delete myfunc; return;}
GeneralFit myfit(myfunc);
myfit.SetDebug(O.lpg);
myfit.SetData(&mydata);
myfit.SetStopChi2(O.stc2);
{for(int i=0;i<myfunc->NPar();i++) {
  char str[10];
  sprintf(str,"P%d",i);
  myfit.SetParam(i,str,Par(i),Step(i),Min(i),Max(i));
}}
if(O.lp>1) myfit.PrintFit();
double c2r = (double) myfit.Fit();
if(O.lp>0) myfit.PrintFit();
if(c2r>0.) {
  c2r = myfit.GetChi2Red();
  cout<<"C2r_Reduit = "<<c2r<<"  nstep="<<myfit.GetNStep()<<endl;
  Vector ParFit(myfunc->NPar());
  for(int i=0;i<myfunc->NPar();i++) ParFit(i)=myfit.GetParm(i);
} else {
  cout<<"echec Fit, rc = "<<c2r<<"  nstep="<<myfit.GetNStep()<<endl;
}

// Mise a disposition des resultats
if(c2r>=0. && myfunc && (O.okres>0||O.okfun>0)) {
  string nomres = nom + "res";
  string nomfun = nom + "fun";
  if(v) {
    if(O.okres) {Vector* ob = v->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {Vector* ob = v->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  } else if(h) {
    if(O.okres) {Histo* ob = h->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {Histo* ob = h->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  } else if(m) {
    if(O.okres) {Matrix* ob = m->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {Matrix* ob = m->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  } else if(h2) {
    if(O.okres) {Histo2D* ob = h2->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {Histo2D* ob = h2->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  } else if(im) {
    if(O.okres) {RzImage* ob = im->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {RzImage* ob = im->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  } else if(g) {
    if(O.okres) {GeneralFitData* ob = g->FitResidus(myfit);  if(ob) mOmg->AddObj(ob,nomres);}
    if(O.okfun) {GeneralFitData* ob = g->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}
  }
}

// Nettoyage
if(myfunc) delete myfunc;
return;
}


/* --Methode-- */
void Services2NObjMgr::ComputeExpressions(NObjMgrAdapter* obja, string& expx, 
                          string& expy, string& expz, string& expt, string& expcut, 
                          NTuple* nt, Histo* h1, Histo2D* h2, HProf* hp)
{
if (obja == NULL) return;
NTupleInterface* objnt = obja->GetNTupleInterface();
if (objnt == NULL)  return;
string  vardec = objnt->VarList_C("_zz6qi_");

PlotExprFunc f = LinkExprFunc(vardec, expx, expy, expz, expt, expcut);
if (!f) { 
  cerr << "Services2NObjMgr::::ComputeExpressions() Error Creation PlotExprFunc " <<  endl; 
  return;
  }

double xnt[10];
float  fxnt[10];

int i,k;
for(i=0; i<10; i++) xnt[i] = 0.;

TRY  {
  double* xn;
  for(k=0; k<objnt->NbLines(); k++)    {
    xn = objnt->GetLineD(k);
    if (f(xn, xnt, xnt+1, xnt+2, xnt+3) != 0) {
      if (nt) { 
         for(i=0; i<4; i++) fxnt[i] = xnt[i]; 
         nt->Fill(fxnt);
         }
      if (h1) h1->Add(xnt[0], xnt[3]);
      if (h2) h2->Add(xnt[0], xnt[1], xnt[3]);
      if (hp) hp->Add(xnt[0], xnt[1], xnt[3]);
      }
    }
 }
CATCH(merr) {
  fflush(stdout); 
  cout << endl; 
  cerr << endl;
  string es = PeidaExc(merr);
  cerr << "Services2NObjMgr::ComputeExpressions()  Exception :" << merr << es;
  } ENDTRY;
 

// Fermeture du fichier .so 
CloseDLL();
return;
}


/* --Methode-- */
PlotExprFunc Services2NObjMgr::LinkExprFunc(string& vardec, string& expx, string& expy, 
                                            string& expz, string& expt, string& cut)
{
FILE *fip;
string fname = TmpDir + "expf_pia_dl.c";
string cmd;
int rc;

cmd = "rm -f " + fname;
rc = system(cmd.c_str());
//DBG printf("LinkExprFunc_Do> %s  (Rc=%d)\n", cmd.c_str(), rc); 

if ((fip = fopen(fname.c_str(), "w")) == NULL)   { 
  string sn = fname; 
  cout << "Services2NObjMgr/LinkExprFunc_Erreur: Pb. Ouverture " << sn << endl; 
  return(NULL); 
  }

// constitution du fichier a compiler
fputs("#include <stdlib.h> \n", fip);
fputs("#include <math.h> \n", fip);
fputs("int expf_pia_dl_func(double* _zz6qi_, double* _rx_6q_, double* _ry_6q_, double* _rz_6q_, double* _rt_6q_) \n{\n", fip);
fprintf(fip,"%s \n", vardec.c_str());
fprintf(fip, "if (!(%s)) { *_rx_6q_ = *_ry_6q_ = *_rz_6q_ = *_rt_6q_ = 0.;  return(0); } \n", cut.c_str());
fprintf(fip, "*_rx_6q_ = %s ; \n", expx.c_str());
fprintf(fip, "*_ry_6q_ = %s ; \n", expy.c_str());
fprintf(fip, "*_rz_6q_ = %s ; \n", expz.c_str());
fprintf(fip, "*_rt_6q_ = %s ; \n", expt.c_str());
fputs("return(1); \n} \n", fip);
fclose(fip);
string func = "expf_pia_dl_func";
return((PlotExprFunc)LinkFunctionFromFile(fname, func));
}


/* --Methode-- */
DlFunction Services2NObjMgr::LinkFunctionFromFile(string const & fname, string const & funcname)
{
//  Le link dynamique
CloseDLL();
dynlink = PDynLinkMgr::BuildFromCFile(fname);
if (dynlink == NULL) { 
  cerr << "Services2NObjMgr/LinkFunctionFromFile_Erreur: Erreur creation/Ouverture SO " << endl; 
  return(NULL); 
  }

DlFunction retfunc = dynlink->GetFunction(funcname);
if (retfunc == NULL) { 
  string sn = funcname; 
  cerr << "Services2NObjMgr/LinkExprFunc_Erreur: Erreur linking " << sn << endl; 
  CloseDLL();
  return(NULL); 
  }
else return(retfunc);
}

/* --Methode-- */
void Services2NObjMgr::CloseDLL()
{
if (dynlink) delete dynlink;   dynlink = NULL;
}


/* --Methode-- */
string Services2NObjMgr::FileName2Name(string const & fn)
{

char fsep[2] = {FILESEP, '\0'};
char tsep[2] = {'.', '\0'};
size_t p = fn.find_last_of(fsep);
size_t l = fn.length();
if (p >= l)  p = 0;
else p++;
size_t q = fn.find_first_of(tsep,p);
if (q < p) q = l;
return(fn.substr(p,q-p));
}


typedef vector<string> GraTok;

/* --Methode-- */
int Services2NObjMgr::DecodeDispOption(string& gratt, bool& fgsrgr)
{
int ropt = Disp_Next;
if (!mImgapp)  return(ropt); 

for(int i=0; i<gratt.length(); i++) gratt[i] = tolower(gratt[i]);

if (fgsrgr) mImgapp->SaveGraphicAtt();

if ( (gratt == "def") || (gratt == "default") )  {  // Remise aux valeurs par defaut = non defini
  mImgapp->SetColAtt();
  mImgapp->SetLineAtt();
  mImgapp->SetFontAtt();
  mImgapp->SetMarkerAtt();
  mImgapp->SetColMapId();
  mImgapp->SetZoomAtt();
  mImgapp->SetAxesAtt();
  mImgapp->SetXYLimits();
  mImgapp->UseXYLimits();
  return(ropt);
  }

// On separe en mots separes par des virgules 
gratt = ","+gratt;
size_t p = 0;
size_t q = 0;
size_t l = gratt.length();
string token;

GraTok grt;

while (q < l)  {
  p = gratt.find_first_not_of(" ,",q+1); // au debut d'un token
  if (p>=l) break;
  q = gratt.find_first_of(" ,",p); // la fin du token;
  token = gratt.substr(p,q-p);
  grt.push_back(token);
  }


static GrAttNames::iterator it;

int k; 
bool fgcont = true;
fgsrgr = false;

for(k=0; k<grt.size(); k++) {
//  cout << "--DBG--SetGraphicAttributes() " << grt[k] << endl;

 // Decodage option affichage (win, next, etc 
  fgcont = true;
  if ( (grt[k] == "win") || (grt[k] == "w") )  ropt =  Disp_Win;
  else if ( (grt[k] == "same") || (grt[k] == "s") )  ropt = Disp_Same;
  else  if ( (grt[k] == "stack") || (grt[k] == "st") )  ropt = Disp_Stack;
  else fgcont = false;
  if (fgcont)   continue;

  // Utilisation limites X-Y
  if ( grt[k] == "xylimits" ) { mImgapp->UseXYLimits(true);  fgsrgr = true;  continue; } 
  // Si c'est une couleur
  it = GrAcolors.find(grt[k]);
  if (it != GrAcolors.end())  { mImgapp->SetColAtt((PIColors)((*it).second.a1)); fgsrgr = true; continue; }
  // Si c'est un attribut de lignes
  it = GrAlines.find(grt[k]);
  if (it != GrAlines.end())  { mImgapp->SetLineAtt((PILineAtt)((*it).second.a1));  fgsrgr = true; continue; }
  // Si c'est un attribut de fontes
  it = GrAfonts.find(grt[k]);
  if (it != GrAfonts.end())  { mImgapp->SetFontAtt((PIFontSize)((*it).second.a2), (PIFontAtt)((*it).second.a1) ); 
                               fgsrgr = true;  continue; }
  // Si c'est un attribut de markers
  it = GrAmarkers.find(grt[k]);
  if (it != GrAmarkers.end())  { mImgapp->SetMarkerAtt((*it).second.a2, (PIMarker)((*it).second.a1) ); 
                                 fgsrgr = true;  continue; }
  // Si c'est un colormap
  it = GrAcmap.find(grt[k]);
  if (it != GrAcmap.end())  { mImgapp->SetColMapId( (CMapId)((*it).second.a1) ); fgsrgr = true; continue; }  
  // Si c'est un facteur de zoom
  it = GrAzoom.find(grt[k]);
  if (it != GrAzoom.end())  { mImgapp->SetZoomAtt( (*it).second.a1 );  fgsrgr = true; continue; }  
  // Si c'est un attribut d'axe
  it = GrAaxes.find(grt[k]);
  if (it != GrAaxes.end())  { mImgapp->SetAxesAtt( (*it).second.a1 );  fgsrgr = true; continue; }  
  
  }  

return(ropt);
}



// !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
// Initialisation des chaines de caracteres designant les attributs graphiques 

/* --Methode-- */
void Services2NObjMgr::InitGrAttNames()
{
gratt_item gi;
// Les couleurs
gi.a2 = 0;
gi.a1 = PI_NotDefColor;
GrAcolors["defcol"] = gi;
gi.a1 = PI_Black;
GrAcolors["black"] = gi;
gi.a1 = PI_White;
GrAcolors["white"] = gi;
gi.a1 = PI_Grey;
GrAcolors["grey"] = gi;
gi.a1 = PI_Red;
GrAcolors["red"] = gi;
gi.a1 = PI_Blue;
GrAcolors["blue"] = gi;
gi.a1 = PI_Green;
GrAcolors["green"] = gi;
gi.a1 = PI_Yellow;
GrAcolors["yellow"] = gi;
gi.a1 = PI_Magenta;
GrAcolors["magenta"] = gi;

gi.a1 = PI_Cyan;
GrAcolors["cyan"] = gi;
gi.a1 = PI_Turquoise;
GrAcolors["turquoise"] = gi;
gi.a1 = PI_NavyBlue;
GrAcolors["navyblue"] = gi;
gi.a1 = PI_Orange;
GrAcolors["orange"] = gi;
gi.a1 = PI_SiennaRed;
GrAcolors["siennared"] = gi;
gi.a1 = PI_Purple;
GrAcolors["purple"] = gi;
gi.a1 = PI_LimeGreen;
GrAcolors["limegreen"] = gi;
gi.a1 = PI_Gold;
GrAcolors["gold"] = gi;

// Les attributs de lignes
gi.a2 = 0;
gi.a1 = PI_NotDefLineAtt;
GrAlines["defline"] = gi;
gi.a1 = PI_NormalLine;
GrAlines["normalline"] = gi;
gi.a1 = PI_ThinLine;
GrAlines["thinline"] = gi;
gi.a1 = PI_ThickLine;
GrAlines["thickline"] = gi;
gi.a1 = PI_DashedLine;
GrAlines["dashedline"] = gi;
gi.a1 = PI_ThinDashedLine;
GrAlines["thindashedline"] = gi;
gi.a1 = PI_ThickDashedLine;
GrAlines["thickdashedline"] = gi;
gi.a1 = PI_DottedLine;
GrAlines["dottedline"] = gi;
gi.a1 = PI_ThinDottedLine;
GrAlines["thindottedline"] = gi;
gi.a1 = PI_ThickDottedLine;
GrAlines["thickdottedline"] = gi;

// Les fontes 
gi.a2 = PI_NotDefFontSize;
gi.a1 = PI_NotDefFontAtt;
GrAlines["deffont"] = gi;

gi.a2 = PI_NormalSizeFont;
gi.a1 = PI_RomanFont;
GrAlines["normalfont"] = gi;
gi.a1 = PI_BoldFont;
GrAlines["boldfont"] = gi;
gi.a1 = PI_ItalicFont;
GrAlines["italicfont"] = gi;
gi.a2 = PI_SmallSizeFont;
gi.a1 = PI_RomanFont;
GrAlines["smallfont"] = gi;
gi.a1 = PI_BoldFont;
GrAlines["smallboldfont"] = gi;
gi.a1 = PI_ItalicFont;
GrAlines["smallitalicfont"] = gi;
gi.a2 = PI_BigSizeFont;
gi.a1 = PI_RomanFont;
GrAlines["bigfont"] = gi;
gi.a1 = PI_BoldFont;
GrAlines["bigboldfont"] = gi;
gi.a1 = PI_ItalicFont;
GrAlines["bigitalicfont"] = gi;
gi.a2 = PI_HugeSizeFont;
gi.a1 = PI_RomanFont;
GrAlines["hugefont"] = gi;
gi.a1 = PI_BoldFont;
GrAlines["hugeboldfont"] = gi;
gi.a1 = PI_ItalicFont;
GrAlines["hugeitalicfont"] = gi;


// Les markers 
const char* mrkn[11] = { "dotmarker", "plusmarker", "crossmarker", 
                         "circlemarker", "fcirclemarker", "boxmarker", "fboxmarker", 
                         "trianglemarker", "ftrianglemarker", "starmarker", "fstarmarker"};
PIMarker mrk[11] = { PI_DotMarker, PI_PlusMarker, PI_CrossMarker, 
                     PI_CircleMarker, PI_FCircleMarker, PI_BoxMarker, PI_FBoxMarker, 
                     PI_TriangleMarker, PI_FTriangleMarker, PI_StarMarker, PI_FStarMarker};

gi.a2 = 0;
gi.a1 = PI_NotDefMarker;
GrAmarkers["defmarker"] = gi;
 
for(int j=0; j<11; j++) {
  string smrk;
  char buff[16];
  for(int m=1; m<10; m+=2) {
    sprintf(buff,"%d",m);
    smrk = (string)mrkn[j] + (string)buff;
    gi.a1 = mrk[j];  gi.a2 = m;
    GrAmarkers[smrk] = gi;
    }
  }

// Les tables de couleurs 
gi.a2 = 0;
gi.a1 = CMAP_OTHER;
GrAcmap["defcmap"] = gi;
gi.a1 = CMAP_GREY32;
GrAcmap["grey32"] = gi;
gi.a1 = CMAP_GREYINV32;
GrAcmap["greyinv32"] = gi;
gi.a1 = CMAP_COLRJ32;
GrAcmap["colrj32"] = gi;
gi.a1 = CMAP_COLBR32;
GrAcmap["colbr32"] = gi;
gi.a1 = CMAP_GREY128;
GrAcmap["grey128"] = gi;
gi.a1 = CMAP_GREYINV128;
GrAcmap["greyinv128"] = gi;
gi.a1 = CMAP_COLRJ128;
GrAcmap["colrj128"] = gi;
gi.a1 = CMAP_COLBR128;
GrAcmap["colbr128"] = gi;

// La valeur de zoom
gi.a2 = 0;
gi.a1 = 0;
GrAzoom["defzoom"] = gi;
gi.a1 = 1;
GrAzoom["zoomx1"] = gi;
gi.a1 = 2;
GrAzoom["zoomx2"] = gi;
gi.a1 = 3;
GrAzoom["zoomx3"] = gi;
gi.a1 = 4;
GrAzoom["zoomx4"] = gi;
gi.a1 = 5;
GrAzoom["zoomx5"] = gi;
gi.a1 = -2;
GrAzoom["zoom/2"] = gi;
gi.a1 = -3;
GrAzoom["zoom/3"] = gi;
gi.a1 = -4;
GrAzoom["zoom/4"] = gi;
gi.a1 = -5;
GrAzoom["zoom/5"] = gi;

// Attributs d'axes
gi.a2 = 0;
gi.a1 = (int)(kBoxAxes | kExtTicks | kLabels);
GrAaxes["stdaxes"] = gi;
GrAaxes["defaxes"] = gi;
GrAaxes["boxaxes"] = gi;
gi.a1 = (int)kAxesDflt;
GrAaxes["simpleaxes"] = gi;
gi.a1 = (int)(kBoxAxes | kExtTicks | kLabels | kGridOn);
GrAaxes["boxaxesgrid"] = gi;

gi.a1 = (int)(kBoxAxes | kTicks | kLabels | kMinTicks | kMajTicks);
GrAaxes["fineaxes"] = gi;
 gi.a1 = (int)(kBoxAxes | kTicks | kLabels | kMinTicks | kMajTicks | kGridOn);
GrAaxes["grid"] = gi;
GrAaxes["fineaxesgrid"] = gi;
                
}


// SANS_EVOLPLANCK     Attention !
#include "pclassids.h"

/* --Methode-- */
char* Services2NObjMgr::PClassIdToClassName(int cid)
{
switch (cid) {
  case ClassId_Poly1 : 
    return("Poly1");
  case ClassId_Poly2 : 
    return("Poly2");
  case ClassId_Matrix : 
    return("Matrix");
  case ClassId_Vector : 
    return("Vector");

  case ClassId_DVList : 
    return("DVList");

  case ClassId_Histo1D : 
    return("Histo1D");
  case ClassId_Histo2D : 
    return("Histo2D");
  case ClassId_HProf : 
    return("HProf");
  case ClassId_NTuple : 
    return("NTuple");
  case ClassId_GeneralFitData : 
    return("GeneralFitData");

  case ClassId_Image : 
    return("RzImage");
  case ClassId_Image + kuint_1 : 
    return("ImageU1");
  case ClassId_Image + kint_1 : 
    return("ImageI1");
  case ClassId_Image + kuint_2 : 
    return("ImageU2");
  case ClassId_Image + kint_2 : 
    return("ImageI2");
  case ClassId_Image + kuint_4 : 
    return("ImageU4");
  case ClassId_Image + kint_4 : 
    return("ImageI4");
  case ClassId_Image + kr_4 : 
    return("ImageR4");
  case ClassId_Image + kr_8 : 
    return("ImageR8");

  case ClassId_ZFidu : 
    return("ZFidu");

  case ClassId_StarList : 
    return("StarList");
  case ClassId_Transfo : 
    return("Transfo");
  case ClassId_PSF : 
    return("PSF");


// - Ajout objet PPF
  default:
    return("AnyDataObj");
  }
}