00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020 #include "EvtGenBase/EvtPatches.hh"
00021
00022 #include <iostream>
00023 #include <fstream>
00024 #include <stdlib.h>
00025 #include <ctype.h>
00026 #include <string.h>
00027 #include "EvtGenBase/EvtOrthogVector.hh"
00028 using std::fstream;
00029
00030 EvtOrthogVector::EvtOrthogVector(int n, std::vector<double> *vectors){
00031
00032 _dimen=n;
00033 _holder.resize(n);
00034
00035 std::vector<int> temp;
00036
00037 int i;
00038 for (i=0;i<n;i++) {
00039 _orthogVector.push_back(0.);
00040 temp.push_back(i);
00041 }
00042
00043 findOrthog(_dimen,temp, vectors);
00044
00045 }
00046
00047 EvtOrthogVector::~EvtOrthogVector(){
00048 }
00049
00050 void EvtOrthogVector::findOrthog(int dim, std::vector<int> invect,
00051 std::vector<double> *vectors) {
00052
00053
00054 if ( dim==2 ) {
00055 _holder[0]=invect[0];
00056 _holder[1]=invect[1];
00057 int sign=findEvenOddSwaps();
00058 {
00059 double addition=1;
00060 int i;
00061 for (i=1; i<_dimen; i++){
00062 addition*=vectors[i-1][_holder[i]];
00063 }
00064 addition*=sign;
00065 _orthogVector[_holder[0]]+=addition;
00066 }
00067
00068 _holder[0]=invect[1];
00069 _holder[1]=invect[0];
00070
00071 {
00072 double addition=1;
00073 int i;
00074 for (i=1; i<_dimen; i++){
00075 addition*=vectors[i-1][_holder[i]];
00076 }
00077 addition*=sign;
00078 _orthogVector[_holder[0]]-=addition;
00079 }
00080
00081 return;
00082 }
00083 else{
00084 std::vector<int> temp((2*dim));
00085
00086 int i;
00087 for (i=0; i<dim; i++) temp[i]=invect[i];
00088 for (i=0; i<dim; i++) temp[i+dim]=invect[i];
00089
00090 for (i=0; i<dim; i++) {
00091 _holder[dim-1]=temp[dim-1+i];
00092 std::vector<int> tempDim((dim-1));
00093
00094 int j;
00095 for (j=0; j<(dim-1); j++) tempDim[j]=temp[j+i];
00096 findOrthog(dim-1, tempDim, vectors);
00097 }
00098 }
00099
00100 return;
00101 }
00102
00103 int EvtOrthogVector::findEvenOddSwaps() {
00104
00105 std::vector<int> temp(_dimen);
00106
00107 int i,j,nSwap;
00108 for (i=0; i<_dimen; i++) temp[i]=_holder[i];
00109
00110 nSwap=0;
00111 for (i=0; i<(_dimen-1); i++) {
00112 for (j=i+1; j<_dimen; j++) {
00113
00114 if ( temp[i]>temp[j] ) {
00115 int duh=temp[j];
00116 temp[j]=temp[i];
00117 temp[i]=duh;
00118 nSwap+=1;
00119 }
00120 }
00121 }
00122 nSwap-= (nSwap/2)*2;
00123
00124 if ( nSwap ) return -1;
00125
00126 return 1;
00127
00128 }