//============================================================// COOOL           version 1.1           ---     Nov,  1995//   Center for Wave Phenomena, Colorado School of Mines//============================================================////   This code is part of a preliminary release of COOOL (CWP// Object-Oriented Optimization Library) and associated class // libraries. //// The COOOL library is a free software. You can do anything you want// with it including make a fortune.  However, neither the authors,// the Center for Wave Phenomena, nor anyone else you can think of// makes any guarantees about anything in this package or any aspect// of its functionality.//// Since you've got the source code you can also modify the// library to suit your own purposes. We would appreciate it // if the headers that identify the authors are kept in the // source code.////==========================================================================// author:  H. Lydia Deng, 08/14/95// Implementation of the Simulaneous Iterative Reconstruction Technique (SIRT)//==========================================================================#include <SIRT.hh>SIRT::~SIRT(){;}SIRT::SIRT(int n, LinearForward* p, Vector<double>* data, 			int it, double eps) : QuadraticOptima(n, p, data){    	iterMax 	= 	it;    	tol 		= 	eps;	m		=	p->modelSize();	n		=	p->dataSize();}SIRT::SIRT(int n, LinearForward* p, Vector<double>* data, 			int it, double eps, int verb) : QuadraticOptima(n, p, data, verb){    	iterMax 	= 	it;    	tol 		= 	eps;	m		=	p->modelSize();	n		=	p->dataSize();}Model<double> SIRT::upDating(int irow, const Model<double>& model){   Model<double> new_model(model);   Vector<double> a(m), data(n), b(m);   MisFitFcn* f = (MisFitFcn*) fp;   double scale, sigma=0;   int i;      data = f->getRefData();   b = 0.;   for (i=0; i<m; i++)   {      a = f->getOneOperator(i);      scale = data[i] - model.modParam()*a;      b += scale*a;      sigma += a.norm2S();   }   sigma = 1./sigma;   new_model = model.update(1,sigma,b);      return new_model;}Model<long> SIRT::upDating(int irow, const Model<long>& model){    Model<double> temp(model);    temp = upDating(irow, temp);    Model<long> new_model(temp);    return new_model;}Model<double> SIRT::optimizer(Model<double>& init_model){   Model<double> model(init_model);   double value = fp->performance(model);   if (isVerbose) cerr << "the current residue: "<< value << endl;   QuadraticOptima::appendResidue(value);   int iter = 0;   while(iter < iterMax)   {      int k = iter%m;      model = upDating(k, model);      iter ++;            if (k == 0)       {	 value = fp->performance(model);	 if (isVerbose) cerr << "the current residue: "<< value << endl;	 QuadraticOptima::appendResidue(value);      }          }    return(model);}Model<long> SIRT::optimizer(Model<long>& model){    Model<double> temp(model);    temp = optimizer(temp);    Model<long> optm(temp);    return optm;}