/*  quick.cc

    Mark Woolrich, FMRIB Image Analysis Group

    Copyright (C) 1999-2000 University of Oxford  */

/*  CCOPYRIGHT  */

#include <iostream>
#include <fstream>
#include <iomanip>
#include <sstream>
#define WANT_STREAM
#define WANT_MATH

#include "newmatap.h"
#include "newmatio.h"
#include <string>
#include <math.h>
#include "glmrand.h"
#include "miscmaths/volume.h"
#include "utils/log.h"
#include "histogram.h"
#include "miscmaths/t2z.h"
#include "miscmaths/f2z.h"
#include "miscmaths/miscmaths.h"
#include "libprob.h"
#include <time.h>

using namespace Utilities;
using namespace NEWMAT;
using namespace FILM;
using namespace MISCMATHS;

int main(int argc, char *argv[])
{
  try{
    rand();

//     Log::getInstance().setDir(".");

//     //    cerr << "z = " << T2z::getInstance().convert(9.02,958) << endl;
//     int X = 20;
//     int Y = 20;
//     int Z = 1;
//     int T = 200;
    
//     int N = X*Y*Z*T;
//     //Matrix mat(N,N);
//     SymmetricBandMatrix mat(N,X*Y*Z);
//     mat = 0;

//     cerr << "X*Y*Z*T=" << N << endl;
//     cerr << "X*Y*Z=" << X*Y*Z << endl;
//     int mba = 0;

//     for(float j=0;j<=0.8;j+=0.4)
//       {
// 	for(float i=-0.8;i<=0.8;i+=0.8)
// 	  {
	
//       for(int t = 1; t <= T; t++)
//       for(int z = 1; z <= Z; z++)
// 	for(int y = 1; y <= Y; y++)
// 	  for(int x = 1; x <= X; x++)
// 	    {
// 	      int a = x+(y-1)*X+(z-1)*X*Y+(t-1)*Y*X*Z;

// 	      for(int t2 = Max(t-1,1); t2 <= Min(t+1,T); t2++)
// 		for(int z2 = Max(z-1,1); z2 <= Min(z+1,Z); z2++)
// 		  for(int y2 = Max(y-1,1); y2 <= Min(y+1,Y); y2++)
// 		    for(int x2 = Max(x-1,1); x2 <= Min(x+1,X); x2++)
// 		      {		    
// 			int b = x2+(y2-1)*X+(z2-1)*X*Y+(t2-1)*Y*X*Z;						
// 			if(b-a > mba)
// 			  mba = b-a;

// // 			cerr <<"t="<<t;
// // 			cerr <<"z="<<z;
// // 			cerr <<"y="<<z;
// // 			cerr <<"x="<<x;
// // 			cerr <<"t2="<<t2;
// // 			cerr <<"z2="<<z2;
// // 			cerr <<"y2="<<z2;
// // 			cerr <<"x2="<<x2<<endl;

// 			if (x==x2 && y==y2 && z==z2 && t==t2)
// 			  {			    
// 			    mat(a,b) = 1;
// 			  }
// 			else
// 			  {
// 			    if(t==t2)
// 			      mat(a,b) = -j;
// 			    else if(x==x2 && y==y2 && z==z2)
// 			      mat(a,b) = -i;
// 			  }			
// 		      }
// 	    }
//       //cerr << "mba=" << mba <<endl;
//       //write_ascii_matrix(mat,"mat");
      
//       clock_t start = clock();
	        
//       //cerr << "det(mat)=" << mat.Determinant() << endl;
      
//       cerr << "logdet(mat)=" << mat.LogDeterminant().LogValue() << endl;
    //	  }
    //      }
      

    //      cerr << "time taken=" << (clock()-start)/float(CLOCKS_PER_SEC) << endl;
  }
  catch(Exception p_excp) 
    {
      cerr << p_excp.what() << endl;
    }
  return 0;
}