source: Sophya/trunk/SophyaLib/SkyT/specrespvector.cc@ 668

//--------------------------------------------------------------------------
// File and Version Information:
//      $Id: specrespvector.cc,v 1.4 1999-11-29 14:16:09 ansari Exp $
//
// Description:
//      Aim of the class: To give the energy density
//                        The unity used here is W/m^2/Hz/sr
//
// History (add to end):
//      Sophie   Oct, 1999  - creation
//
//------------------------------------------------------------------------

//---------------
// C++ Headers --
//---------------
#include "machdefs.h"
#include <iostream.h>
#include <math.h>
// #include <typeinfo>

#include "specrespvector.h"
#include "pexceptions.h"
//----------------
// Constructor --
//----------------
SpecRespVec::SpecRespVec()
        : SpectralResponse()
{
}

SpecRespVec::SpecRespVec(Vector const & nu, Vector const & fdenu, double numin, double numax)
        : SpectralResponse(numin, numax)
{
  if(nu.NElts() != fdenu.NElts())  
    throw SzMismatchError("SpecRespVec::SpecRespVec() - Non equal vector sizes");
  _vecOfNu = nu;
  _vecOfFDeNu = fdenu;
  _numin = nu(0);
  _numax = nu(nu.NElts()-1);
  _size = nu.NElts();
  if(_vecOfNu.NElts() != _vecOfFDeNu.NElts()) cout << "vectors are not compatible" << exit;
}


//--------------
// Destructor --
//--------------
SpecRespVec::~SpecRespVec()
{
}

//              ---------------------------
//              --  Function Definitions --
//              ---------------------------



 
double 
SpecRespVec::transmission(double nu)  const 
{
  if ( (nu < _numin) || (nu > _numax) ) return(0.);
  double value = 0.;
  int sizeVecOfNu = _vecOfNu.NElts();
  if(nu <=  _vecOfNu(0)) return _vecOfFDeNu(0);
  if(nu >=  _vecOfNu(sizeVecOfNu-1)) return _vecOfFDeNu(sizeVecOfNu-1);

  for (int i=1; i<sizeVecOfNu; i++)
    {
      if(nu < _vecOfNu(i))
        {
          double up = _vecOfFDeNu(i) ;
          double down = _vecOfFDeNu(i-1);
          double xmin = _vecOfNu(i-1);
          double xmax = _vecOfNu(i);
          double a = ((up-down)/(xmax-xmin));
          value = a*nu+(up-a*xmax);
          return value;
        }
    }
  return value;
}



void
SpecRespVec::Print(ostream& os) const 
{

  //   os << "SpecRespVec::Print (" << typeid(*this).name() 
  //<< ") - Fmin,Fmax= " << minFreq() << "," << maxFreq() << endl;
  os << "SpecRespVec ::Print - Fmin,Fmax= " << minFreq() << "," << maxFreq() << endl;
  os << "MeanFreq= " << meanFreq() << "  Transmission= " << transmission(meanFreq()) << endl;
  os << "PeakFreq= " << peakFreq() << "  Transmission= " << transmission(peakFreq()) << endl;

}


void
ObjFileIO<SpecRespVec>::WriteSelf(POutPersist& s) const
{
  if(dobj == NULL)
    {
      cout << " ObjFileIO<SpecRespVec>::WriteSelf:: dobj= null " << endl;
      return;
    }

  int_4 version, nothing;
  version = 1;
  nothing = 0;   // Reserved for future use 
  s.PutI4(version);
  s.PutI4(nothing);

  s.PutR8(dobj->minFreq());
  s.PutR8(dobj->maxFreq());

  // TVector<T> has Persistence Manager 
  s <<  dobj->getNuVec();
  s <<  dobj->getTNuVec();
}

void 
ObjFileIO<SpecRespVec>::ReadSelf(PInPersist& s) 
{
  int_4 version, nothing;
  version = 1;
  nothing = 0;   // Reserved for future use 
  s.GetI4(version);
  s.GetI4(nothing);

  if(dobj == NULL)
    {
      Vector v1(10);
      Vector v2(10);
      dobj= new SpecRespVec();
      ownobj= true;
    }
  
  r_8 minf, maxf;
  s.GetR8(minf);
  s.GetR8(maxf);
  dobj->setMinMaxFreq(minf, maxf);
  // TVector<T> has Persistence Manager 
  s >>  dobj->getNuVec();
  s >>  dobj->getTNuVec();
}


#ifdef __CXX_PRAGMA_TEMPLATES__
#pragma define_template ObjFileIO<SpecRespVec>
#endif

#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
template class ObjFileIO<SpecRespVec>;
#endif