//
//  XSPEC12  November 2003
//
//

#include <funcWrappers.h>
#include <functionMap.h>
#include <xsTypes.h>
#include <cfortran.h>

void cppModelWrapper(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr, 
        int nPar, void (*cppFunc)(const RealArray&, const RealArray&,
        int, RealArray&, RealArray&, const string&));

void fcppModelWrapper(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError,  
        int nPar, void (*cppFunc)(const RealArray&, const RealArray&,
	int, RealArray&, RealArray&, const string&));

void f_powerLaw(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_zpowerLaw(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_cutoffPowerLaw(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_cemVMekal(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_kdblur2(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux); 

void f_kerrbb(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_kerrdisk(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_pileup(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_xneq(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux); 

void f_xsbexriv(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux); 

void f_xslaor(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_xsneq(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux); 

void f_xspexriv(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux); 

void f_cflux(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_PartialCovering(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_simpl(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

void f_superExpCutoff(const float* energy, int nFlux, const float* params, 
	     int spectrumNumber, float* flux, float* fluxError); 

FCALLSCSUB6(f_powerLaw,XSPWLW,xspwlw,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_zpowerLaw,XSZPLW,xszplw,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_cutoffPowerLaw,XSPLCO,xsplco,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_cemVMekal,CEVMKL,cevmkl,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB5(f_kdblur2,KDBLUR2,kdblur2,FLOATV,INT,FLOATV,INT,FLOATV) 
FCALLSCSUB6(f_kerrbb,KERRBB,kerrbb,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_kerrdisk,KERRDISK,kerrdisk,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_pileup,PILEUP,pileup,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB5(f_xneq,XNEQ,xneq,FLOATV,INT,FLOATV,INT,FLOATV)
FCALLSCSUB5(f_xsbexriv,XSBEXRIV,xsbexriv,FLOATV,INT,FLOATV,INT,FLOATV)
FCALLSCSUB6(f_xslaor,XSLAOR,xslaor,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB5(f_xsneq,XSNEQ,xsneq,FLOATV,INT,FLOATV,INT,FLOATV)
FCALLSCSUB5(f_xspexriv,XSPEXRIV,xspexriv,FLOATV,INT,FLOATV,INT,FLOATV)
FCALLSCSUB6(f_cflux,CFLUX,cflux,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_PartialCovering,PRTCOV,prtcov,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_simpl,SIMPL,simpl,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 
FCALLSCSUB6(f_superExpCutoff,SPEXPCUT,spexpcut,FLOATV,INT,FLOATV,INT,FLOATV,FLOATV) 

FCALLSCSUB7(C_powerLaw,CXSPWLW,cxspwlw,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_zpowerLaw,CXSZPLW,cxszplw,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_cutoffPowerLaw,CXSPLCO,cxsplco,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_cemVMekal,CCEVMKL,ccevmkl,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_kdblur2,CKDBLU2,ckdblu2,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_kerrbb,CKERRBB,ckerrbb,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_kerrdisk,CKERRDIS,ckerrdis,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_pileup,CPILEUP,cpileup,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_xneq,CXNEQ,cxneq,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_xsbexriv,CXSBEXRI,cxsbexri,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_xslaor,CXSLAOR,cxslaor,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_xsneq,CXSNEQ,cxsneq,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_xspexriv,CXSPEXRI,cxspexri,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_cflux,CCFLUX,ccflux,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_PartialCovering,CPRTCOV,cprtcov,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_simpl,CSIMPL,csimpl,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 
FCALLSCSUB7(C_superExpCutoff,CSPEXPCU,cspexpcu,DOUBLEV,INT,DOUBLEV,INT,DOUBLEV,DOUBLEV,STRING) 



void C_powerLaw(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 1;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, powerLaw);
}

void C_zpowerLaw(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 2;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, zpowerLaw);
}

void C_cutoffPowerLaw(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 2;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, cutoffPowerLaw);
}

void C_cemVMekal(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 19;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, cemVMekal);
}

void C_kdblur2(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 6;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, kdblur2);
}

void C_kerrbb(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 9;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, kerrbb);
}

void C_kerrdisk(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 7;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, kerrdisk);
}

void C_pileup(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 6;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, pileup);
}

void C_xneq(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 4;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, xneq);
}

void C_xsbexriv(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 11;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, xsbexriv);
}

void C_xslaor(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 5;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, xslaor);
}

void C_xsneq(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 16;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, xsneq);
}

void C_xspexriv(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 9;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, xspexriv);
}

void C_cflux(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 3;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, cflux);
}

void C_PartialCovering(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 1;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, PartialCovering);
}

void C_simpl(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 3;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, simpl);
}

void C_superExpCutoff(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr) 
{
   const size_t nPar = 2;
   cppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        initStr, nPar, superExpCutoff);
}


void f_powerLaw(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 1;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, powerLaw);
}

void f_zpowerLaw(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 2;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, zpowerLaw);
}

void f_cutoffPowerLaw(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 2;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, cutoffPowerLaw);
}

void f_cemVMekal(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 19;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, cemVMekal);
}

void f_kdblur2(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux) 
{
   const size_t nPar = 6;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, 0,
        nPar, kdblur2);
}

void f_kerrbb(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 9;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, kerrbb);
}

void f_kerrdisk(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 7;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, kerrdisk);
}

void f_pileup(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 6;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, pileup);
}

void f_xneq(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux) 
{
   const size_t nPar = 4;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, 0,
        nPar, xneq);
}

void f_xsbexriv(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux) 
{
   const size_t nPar = 11;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, 0,
        nPar, xsbexriv);
}

void f_xslaor(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 5;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, xslaor);
}

void f_xsneq(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux) 
{
   const size_t nPar = 16;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, 0,
        nPar, xsneq);
}

void f_xspexriv(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux) 
{
   const size_t nPar = 9;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, 0,
        nPar, xspexriv);
}

void f_cflux(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 3;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, cflux);
}

void f_PartialCovering(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 1;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, PartialCovering);
}

void f_simpl(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 3;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, simpl);
}

void f_superExpCutoff(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError) 
{
   const size_t nPar = 2;
   fcppModelWrapper(energy, nFlux, params, spectrumNumber, flux, fluxError,
        nPar, superExpCutoff);
}


void cppModelWrapper(const double* energy, int nFlux, const double* params, 
        int spectrumNumber, double* flux, double* fluxError, const char* initStr, 
        int nPar, void (*cppFunc)(const RealArray&, const RealArray&,
        int, RealArray&, RealArray&, const string&))
{
   // Assumes energy points to arrays of size nFlux+1, flux and fluxError
   // point to arrays of size nFlux (though they need not be initialized),
   // and params points to an array of size nPar.
   RealArray energy_C(energy, (size_t)nFlux+1);
   RealArray params_C(params, nPar);
   RealArray flux_C(flux, (size_t)nFlux);
   RealArray fluxError_C(fluxError, (size_t)nFlux);
   string cppStr;
   if(initStr && strlen(initStr))
      cppStr = initStr;
   (*cppFunc)(energy_C, params_C, spectrumNumber, flux_C, fluxError_C, cppStr);
   for (int i=0; i<nFlux; ++i)
   {
      flux[i] = flux_C[i];
   }
   if (fluxError_C.size())
   {
      for (int i=0; i<nFlux; ++i)
      {
         fluxError[i] = fluxError_C[i];
      }
   }
} 

void fcppModelWrapper(const float* energy, int nFlux, const float* params, 
        int spectrumNumber, float* flux, float* fluxError, 
        int nPar, void (*cppFunc)(const RealArray&, const RealArray&,
        int, RealArray&, RealArray&, const string&))
{
   // Assumes energy points to arrays of size nFlux+1, flux and fluxError
   // point to arrays of size nFlux (though they need not be initialized),
   // and params points to an array of size nPar.
   RealArray energy_C(0.0, (size_t)nFlux+1);
   RealArray params_C(0.0, nPar);
   RealArray flux_C(0.0, (size_t)nFlux);
   RealArray fluxError_C(0.0, (size_t)nFlux);
   string cppStr;

   for (int i=0; i<nFlux+1; ++i)
   {
      energy_C[i] = (double)energy[i];
   }
   for (int i=0; i<nPar; ++i)
   {
      params_C[i] = (double)params[i];
   }
   for (int i=0; i<nFlux; ++i)
   {
      flux_C[i] = (double)flux[i];
   }
   if (fluxError)
   {
      for (int i=0; i<nFlux; ++i)
      {
         fluxError_C[i] = (double)fluxError[i];
      }
   }


   (*cppFunc)(energy_C, params_C, spectrumNumber, flux_C, fluxError_C, cppStr);

   for (int i=0; i<nFlux; ++i)
   {
      flux[i] = (float)flux_C[i];
   }
   if (fluxError_C.size())
   {
      for (int i=0; i<nFlux; ++i)
      {
         fluxError[i] = (float)fluxError_C[i];
      }
   }
}