GPLIB++: disttest.cpp Source File

00001 #include "MTTensor.h"
00002 #include <fstream>
00003 #include <iostream>
00004 #include <complex>
00005 #include "UniformRNG.h"
00006 
00007 using namespace std;
00008 int main()
00009 {
00010         const complex<double> xx(1.405,2.23);
00011         const complex<double> xy(5.33,2.5);
00012         const complex<double> yx(-7.45,-4.23);
00013         const complex<double> yy(-1.45,-3.29);
00014         const MTTensor TestData(xx,xy,yx,yy);
00015         const double noiselevel = 0.0;
00016         cerr << TestData.GetPhi11() << " " <<  TestData.GetPhi12() << " " <<  TestData.GetPhi21() << " " <<  TestData.GetPhi22() << endl;
00017         cerr << TestData.GetPhiMax() << " " << TestData.GetPhiMin() <<  " " << TestData.GetAlpha_phi() <<  " " << TestData.GetBeta_phi() << endl;
00018         UniformRNG Random;
00019         ofstream zxxfile("zxx.out");
00020         const size_t ntestcases = 10000;
00021         for (size_t i = 0; i < ntestcases; ++i)
00022         {
00023                 double axx = Random.GetNumber() * 100;
00024                 double axy = Random.GetNumber() * 100;
00025                 double ayx = Random.GetNumber() * 100;
00026                 double ayy = Random.GetNumber() * 100;
00027                 complex<double> noisexx(xx.real()+xx.real() *Random.GetNumber() * noiselevel, 
00028                                                 xx.imag() + xx.imag() * Random.GetNumber() * noiselevel);
00029                 complex<double> noisexy(xy.real()+xy.real() *Random.GetNumber() * noiselevel, 
00030                                                 xy.imag() + xy.imag() * Random.GetNumber() * noiselevel);
00031                 complex<double> noiseyx(yx.real()+yx.real() *Random.GetNumber() * noiselevel, 
00032                                                 yx.imag() + yx.imag() * Random.GetNumber() * noiselevel);
00033                 complex<double> noiseyy(yy.real()+yy.real() *Random.GetNumber() * noiselevel, 
00034                                                 yy.imag() + yy.imag() * Random.GetNumber() * noiselevel);
00035                 complex<double> newxx = axx*noisexx+axy*noiseyx;
00036                 complex<double> newxy = axx*noisexy+axy*noiseyy;
00037                 complex<double> newyx = ayx*noisexx+ayy*noiseyx;
00038                 complex<double> newyy = ayx*noisexy+ayy*noiseyy;
00039                 zxxfile << noisexx.real() << " " << noisexx.imag() << endl;
00040                 const MTTensor DistData(newxx,newxy,newyx,newyy);
00041                 //cout << xx << " " << xy << " " << yx << " " << yy << endl;
00042                 //cout << axx << " " << axy << " " << ayx << " " << ayy << endl;
00043                 //cout << newxx << " " << newxy << " " << newyx << " " << newyy << endl;
00044                 cout << (TestData.GetPhiMax() - DistData.GetPhiMax())/TestData.GetPhiMax() << endl;
00045                 
00046         }
00047 }