diff --git a/.AD.cxx.swp b/.AD.cxx.swp
new file mode 100644
index 0000000..bf9330b
Binary files /dev/null and b/.AD.cxx.swp differ
diff --git a/AD.cxx b/AD.cxx
index 9a1d762..162f952 100644
--- a/AD.cxx
+++ b/AD.cxx
@@ -16,780 +16,788 @@ using namespace std;
 int main(int argc,char **argv){
 
 
-		double j1, j2;
+    double j1, j2;
 
 
-		//detector data input
-		double detradius, targetdistance, detthickness = -1;
-		double Energy;
-		//input .csv file
-		double Sigma; 	 
+    //detector data input
+    double detradius, targetdistance, detthickness = -1;
+    double Energy;
+    //input .csv file
+    double Sigma; 	 
 
-		//--------------------------------
+    //--------------------------------
 
-		detradius = 3.;
-		targetdistance = 4.;
-		detthickness = 5.;
+    detradius = 3.;
+    targetdistance = 4.;
+    detthickness = 5.;
 
-		Energy = 1062.;
-		Sigma = 0.; 
+    Energy = 1062.;
+    Sigma = 0.5; 
 
 
-		j1 = 2.; 
-		j2 = 1.;
-		//--------------------------------
+    j1 = 2.; 
+    j2 = 1.;
+    //--------------------------------
 
-		//TEST MODE? y : 1 || n = 0
+    //TEST MODE? y : 1 || n = 0
 
-		int test = 1; 
+    int test = 1; 
 
-		int det_param_token = -1;
-		int gamma_energy_token = -1;
-		int ang_file_token = -1; 
-		int sigma_token = -1; 
+    int det_param_token = -1;
+    int gamma_energy_token = -1;
+    int ang_file_token = -1; 
+    int sigma_token = -1; 
 
-		int j1j2token = -1; 
-		if(test == 1){
+    int j1j2token = -1; 
+    if(test == 1){
 
-				det_param_token = 1;
-				gamma_energy_token = 1;
-				ang_file_token = -1; 
-				sigma_token = 1; 
+        det_param_token = 1;
+        gamma_energy_token = 1;
+        ang_file_token = -1; 
+        sigma_token = 1; 
 
-				j1j2token = 1; 
+        j1j2token = 1; 
 
-		}
-		if( test == 0 ) {
-				printf("Input detector radius, target distance, and detector thickness : ");
-				scanf("%lf,%lf,%lf",&detradius,&targetdistance,&detthickness);
+    }
+    if( test == 0 ) {
+        printf("Input detector radius, target distance, and detector thickness : ");
+        scanf("%lf,%lf,%lf",&detradius,&targetdistance,&detthickness);
 
 
-				if(detradius > 0 && targetdistance > 0 && detthickness > 0){
-						det_param_token = 1;
-				}else{
+        if(detradius > 0 && targetdistance > 0 && detthickness > 0){
+            det_param_token = 1;
+        }else{
 
-						do{
-								if(detradius< 0|| targetdistance < 0|| detthickness < 0){
-										printf("Negative numbers are not allowed!\nRe-enter radius, distance, and thickness : ");
-										scanf("%lf,%lf,%lf", &detradius,&targetdistance,&detthickness);
-								}
-						}while(detradius< 0|| targetdistance < 0|| detthickness < 0);}
+            do{
+                if(detradius< 0|| targetdistance < 0|| detthickness < 0){
+                    printf("Negative numbers are not allowed!\nRe-enter radius, distance, and thickness : ");
+                    scanf("%lf,%lf,%lf", &detradius,&targetdistance,&detthickness);
+                }
+            }while(detradius< 0|| targetdistance < 0|| detthickness < 0);}
 
-				if(detradius > 0 && targetdistance > 0 && detthickness > 0){
-						det_param_token = 1; printf(" --- Detector characteristics loaded --- \n");
-				}
+        if(detradius > 0 && targetdistance > 0 && detthickness > 0){
+            det_param_token = 1; printf(" --- Detector characteristics loaded --- \n");
+        }
 
-				printf("Enter the gamma-ray energy : ");
-				scanf("%lf", &Energy);
-				if(Energy > 0){
-						gamma_energy_token = 1;
-				}else{
-						do{
-								if(Energy <= 0){		
-										printf("Negative numbers are not allowed!\nRe-enter gamma-ray energy : ");
-										scanf("%lf", &Energy);}
+        printf("Enter the gamma-ray energy : ");
+        scanf("%lf", &Energy);
+        if(Energy > 0){
+            gamma_energy_token = 1;
+        }else{
+            do{
+                if(Energy <= 0){		
+                    printf("Negative numbers are not allowed!\nRe-enter gamma-ray energy : ");
+                    scanf("%lf", &Energy);}
 
-						}while(Energy <= 0);}
+            }while(Energy <= 0);}
 
-				if(Energy > 0){
-						gamma_energy_token = 1; printf(" --- Gamma Energy loaded --- \n");
-				}
+        if(Energy > 0){
+            gamma_energy_token = 1; printf(" --- Gamma Energy loaded --- \n");
+        }
 
-				//then calucate QD2 and QD4, and replace the 0 with them. 
-		}
+        //then calucate QD2 and QD4, and replace the 0 with them. 
+    }
 
-		double QD2 =  QK2(Energy, detradius, targetdistance, detthickness); ;
-		//.82;//NEED TO READ FROM ad.txt in the future. 
-		double QD4 =  QK4(Energy, detradius, targetdistance, detthickness); ;
-		//.44;
-		//	cout << "QD2 = " << " " << QD2 << "  QD4 = " << " " << QD4 << "\n";
+    double QD2 =  QK2(Energy, detradius, targetdistance, detthickness); ;
+    //.82;//NEED TO READ FROM ad.txt in the future. 
+    double QD4 =  QK4(Energy, detradius, targetdistance, detthickness); ;
+    //.44;
+    //	cout << "QD2 = " << " " << QD2 << "  QD4 = " << " " << QD4 << "\n";
 
-		//input data file of Angular Data (theta, Yexp, Yerr)
+    //input data file of Angular Data (theta, Yexp, Yerr)
 
-		vector<vector<string> > content;
-		vector<string> row;
-		string line, word; 
+    vector<vector<string> > content;
+    vector<string> row;
+    string line, word; 
 
-		vector<string> adata;
-		vector<string> ydata;
-		vector<string> eydata; 
-		int num = 0;
-		string fname,fname2;
-		cout<<"Enter the file name : ";
-		cin>>fname;
+    vector<string> adata;
+    vector<string> ydata;
+    vector<string> eydata; 
+    int num = 0;
+    string fname,fname2;
+    cout<<"Enter the file name : ";
+    cin>>fname;
 
-		int aaaa = -1;
-		fstream file(fname.c_str());
+    int aaaa = -1;
+    fstream file(fname.c_str());
 
-		if(file.is_open()){
-				cout<< " --- Angular Data File Loaded --- \n";
-				ang_file_token = 1;
-				while(getline(file,line)){
+    if(file.is_open()){
+        cout<< " --- Angular Data File Loaded --- \n";
+        ang_file_token = 1;
+        while(getline(file,line)){
 
-						row.clear();
+            row.clear();
 
-						stringstream str(line);
+            stringstream str(line);
 
-						while(getline(str, word, ','))
+            while(getline(str, word, ','))
 
-								row.push_back(word);
-						content.push_back(row);
-						adata.push_back(row[0]);
-						ydata.push_back(row[1]);
-						eydata.push_back(row[2]);			
+                row.push_back(word);
+            content.push_back(row);
+            adata.push_back(row[0]);
+            ydata.push_back(row[1]);
+            eydata.push_back(row[2]);			
 
-				}
-		}else{
-				do{
-						cout<< "Could not open the file\n";
-						cout<< "Enter the file name : ";
-						cin>>fname2;
+        }
+    }else{
+        do{
+            cout<< "Could not open the file\n";
+            cout<< "Enter the file name : ";
+            cin>>fname2;
 
-						fstream file(fname2.c_str());
+            fstream file(fname2.c_str());
 
-						if(file.is_open()){
-								aaaa = 1;
+            if(file.is_open()){
+                aaaa = 1;
 
-								while(getline(file,line)){
+                while(getline(file,line)){
 
-										row.clear();
+                    row.clear();
 
-										stringstream str(line);
+                    stringstream str(line);
 
-										while(getline(str, word, ','))
+                    while(getline(str, word, ','))
 
-												row.push_back(word);
-										content.push_back(row);
-										adata.push_back(row[0]);
-										ydata.push_back(row[1]);
-										eydata.push_back(row[2]);			
+                        row.push_back(word);
+                    content.push_back(row);
+                    adata.push_back(row[0]);
+                    ydata.push_back(row[1]);
+                    eydata.push_back(row[2]);			
 
-								}
-
-						}else{ aaaa = -1; }
+                }
+
+            }else{ aaaa = -1; }
 
-				}while(aaaa == -1);
+        }while(aaaa == -1);
 
-				cout<< " --- Angular Data File Loaded --- \n";
-				ang_file_token = 1;
+        cout<< " --- Angular Data File Loaded --- \n";
+        ang_file_token = 1;
 
-		}
+    }
 
-		for(int i = 0; i< adata.size();i++){
-				//		cout << adata[i] << "\n";
-				//		cout << ydata[i] << "\n";
-				//		cout << eydata[i] << "\n";
-		}
-
-		// Legendre Polinomial fit using adata, ydata, eydata;
-
-		adata.erase(adata.begin());
-		ydata.erase(ydata.begin());
-		eydata.erase(eydata.begin());
-
-		vector<double> dangle = string_to_double_vector(adata);
-		vector<double> dydata = string_to_double_vector(ydata);
-		vector<double> deydata = string_to_double_vector(eydata);
-
-		vector<double> dangler;
-
-		for(int i = 0; i < dangle.size(); i++){
-				double aa = dangle[i]; 
-
-				dangler.push_back(aa*3.14159/180.);
-				printf("dangle = %lf\n",dangler[i]);
-		}
-		// if you need to scale the y data by sin(theta)	
-
-
-		/*
-		   for(int i=0; i<content.size(); i++){
-		   for(int j=0; j<content[i].size(); j++){
-		   cout<< content[i][j]<<" ";
-		   }
-		   cout<<"\n";
-		   }
-		 */
-		//legendre fitting 
-		// F[xi] = A0(1) + A2(P2(cos(xi))) + A4(P4(cos(xi)))
-		//xi data is angle data in radians. 
-		//yi data is y-intensity data. 
-
-		//to change the number of angles, just mess with the indexing of the read in and arrays. 
-
-
-		//12 constants. this will build our matrix for gaussian elinmination. 
-		//[a1 a2   a3] [A0] = [a4]
-		//[a5 a6   a7] [A2] = [a8]
-		//[a9 a10 a11] [A4] = [a12]
-
-		// Note, the y data must be scaled by sin(theta) of the given angle. 
+    for(int i = 0; i< adata.size();i++){
+        //		cout << adata[i] << "\n";
+        //		cout << ydata[i] << "\n";
+        //		cout << eydata[i] << "\n";
+    }
+
+    // Legendre Polinomial fit using adata, ydata, eydata;
+
+    adata.erase(adata.begin());
+    ydata.erase(ydata.begin());
+    eydata.erase(eydata.begin());
+
+    vector<double> dangle = string_to_double_vector(adata);
+    vector<double> dydata = string_to_double_vector(ydata);
+    vector<double> deydata = string_to_double_vector(eydata);
+
+    vector<double> dangler;
+
+    for(int i = 0; i < dangle.size(); i++){
+        double aa = dangle[i]; 
+
+        dangler.push_back(aa*3.14159/180.);
+        printf("dangle = %lf\n",dangler[i]);
+    }
+    // if you need to scale the y data by sin(theta)	
+
+
+    /*
+       for(int i=0; i<content.size(); i++){
+       for(int j=0; j<content[i].size(); j++){
+       cout<< content[i][j]<<" ";
+       }
+       cout<<"\n";
+       }
+       */
+    //legendre fitting 
+    // F[xi] = A0(1) + A2(P2(cos(xi))) + A4(P4(cos(xi)))
+    //xi data is angle data in radians. 
+    //yi data is y-intensity data. 
+
+    //to change the number of angles, just mess with the indexing of the read in and arrays. 
+
+
+    //12 constants. this will build our matrix for gaussian elinmination. 
+    //[a1 a2   a3] [A0] = [a4]
+    //[a5 a6   a7] [A2] = [a8]
+    //[a9 a10 a11] [A4] = [a12]
+
+    // Note, the y data must be scaled by sin(theta) of the given angle. 
 
-		double a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12 = 0.;
-		/*	
-			for(int i = 0; i< dangler.size(); i++){
-
-		//eq1
-		a1 += 1; //how many data points, not intensional but convientient.  
-		a2 += (1.5 * pow(dangler[i],2) - .5);
-		a3 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
-		a4 += dydata[i]; 
-		//eq2 
-		a5 += (1.5 * pow(dangler[i],2) - .5);
-		a6 += (1.5 * pow(dangler[i],2) - .5)*(1.5 * pow(dangler[i],2) - .5);
-		a7 += (1.5 * pow(dangler[i],2) - .5)*(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
-		a8 += (1.5 * pow(dangler[i],2) - .5)*dydata[i];
-		//eq4
-		a9 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
-		a10 += (1.5 * pow(dangler[i],2) - .5)*(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
-		a11 +=(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. )* (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
-		a12 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. )*dydata[i];
+    double a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12 = 0.;
+    /*	
+        for(int i = 0; i< dangler.size(); i++){
+
+    //eq1
+    a1 += 1; //how many data points, not intensional but convientient.  
+    a2 += (1.5 * pow(dangler[i],2) - .5);
+    a3 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
+    a4 += dydata[i]; 
+    //eq2 
+    a5 += (1.5 * pow(dangler[i],2) - .5);
+    a6 += (1.5 * pow(dangler[i],2) - .5)*(1.5 * pow(dangler[i],2) - .5);
+    a7 += (1.5 * pow(dangler[i],2) - .5)*(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
+    a8 += (1.5 * pow(dangler[i],2) - .5)*dydata[i];
+    //eq4
+    a9 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
+    a10 += (1.5 * pow(dangler[i],2) - .5)*(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
+    a11 +=(35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. )* (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. );
+    a12 += (35./8. * pow(dangler[i],4) - 30./8. * pow(dangler[i],2)  +  3./8. )*dydata[i];
 
 
-		}
+    }
 
-		 */
+*/
 
-		for(int i = 0; i< dangler.size(); i++){
+    for(int i = 0; i< dangler.size(); i++){
 
-				//eq1
-				a1 += 1; //how many data points, not intensional but convientient.  
-				a2 += (1.5 * pow(cos(dangler[i]),2) - .5);
-				a3 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
-				a4 += dydata[i]; 
-				//eq2 
-				a5 += (1.5 * pow(cos(dangler[i]),2) - .5);
-				a6 += (1.5 * pow(cos(dangler[i]),2) - .5)*(1.5 * pow(cos(dangler[i]),2) - .5);
-				a7 += (1.5 * pow(cos(dangler[i]),2) - .5)*(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
-				a8 += (1.5 * pow(cos(dangler[i]),2) - .5)*dydata[i];
-				//eq4
-				a9 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
-				a10 += (1.5 * pow(cos(dangler[i]),2) - .5)*(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
-				a11 +=(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. )* (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
-				a12 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. )*dydata[i];
-
-
-		}
-		/*
-		   printf("a1 = %lf\n",a1);	
-		   printf("a2 = %lf\n",a2);	
-		   printf("a3 = %lf\n",a3);	
-		   printf("a4 = %lf\n",a4);	
-		   printf("a5 = %lf\n",a5);	
-		   printf("a6 = %lf\n",a6);	
-		   printf("a7 = %lf\n",a7);	
-		   printf("a8 = %lf\n",a8);	
-		   printf("a9 = %lf\n",a9);	
-		   printf("a10 = %lf\n",a10);	
-		   printf("a11 = %lf\n",a11);	
-		   printf("a12 = %lf\n",a12);	
-		 */
-		//here number of equations are 3 because its a 3 coeff. fit. 
-		int n = 3; 
-		double mat[3][4] = {{a1,a2,a3,a4},{a5,a6,a7,a8},{a9,a10,a11,a12}};
-		//where A0,A2,A4 are stored
-		double residual[n];
-
-		for(int i = 0; i<n;i++){
-				for( int j = i+1; j<n; j++){
-						if(abs(mat[i][i]) < abs(mat[i][j])){
-								for(int k = 0; k<n+1;k++){
-										//swapping matrix
-										//			mat[i][k] = mat[i][k] + mat[j][k];
-										//			mat[j][k] = mat[i][k] - mat[j][k];
-										//			mat[i][k] = mat[i][k] - mat[j][k];
-										swap(mat[i][k],mat[j][k]);
-								}
-
-						}
-				}
-		}
-		//performing gaussian elimination.
-
-		for(int i = 0; i< n-1; i++) {
-
-				for(int j = i +1; j<n;j++){
-
-						float f = mat[j][i]/mat[i][i];
-
-						for(int k = 0; k<n+1; k++){
-								mat[j][k] = mat[j][k] -f*mat[i][k];
+        //eq1
+        a1 += 1; //how many data points, not intensional but convientient.  
+        a2 += (1.5 * pow(cos(dangler[i]),2) - .5);
+        a3 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
+        a4 += dydata[i]; 
+        //eq2 
+        a5 += (1.5 * pow(cos(dangler[i]),2) - .5);
+        a6 += (1.5 * pow(cos(dangler[i]),2) - .5)*(1.5 * pow(cos(dangler[i]),2) - .5);
+        a7 += (1.5 * pow(cos(dangler[i]),2) - .5)*(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
+        a8 += (1.5 * pow(cos(dangler[i]),2) - .5)*dydata[i];
+        //eq4
+        a9 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
+        a10 += (1.5 * pow(cos(dangler[i]),2) - .5)*(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
+        a11 +=(35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. )* (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. );
+        a12 += (35./8. * pow(cos(dangler[i]),4) - 30./8. * pow(cos(dangler[i]),2)  +  3./8. )*dydata[i];
+
+
+    }
+    /*
+       printf("a1 = %lf\n",a1);	
+       printf("a2 = %lf\n",a2);	
+       printf("a3 = %lf\n",a3);	
+       printf("a4 = %lf\n",a4);	
+       printf("a5 = %lf\n",a5);	
+       printf("a6 = %lf\n",a6);	
+       printf("a7 = %lf\n",a7);	
+       printf("a8 = %lf\n",a8);	
+       printf("a9 = %lf\n",a9);	
+       printf("a10 = %lf\n",a10);	
+       printf("a11 = %lf\n",a11);	
+       printf("a12 = %lf\n",a12);	
+       */
+    //here number of equations are 3 because its a 3 coeff. fit. 
+    int n = 3; 
+    double mat[3][4] = {{a1,a2,a3,a4},{a5,a6,a7,a8},{a9,a10,a11,a12}};
+    //where A0,A2,A4 are stored
+    double residual[n];
+
+    for(int i = 0; i<n;i++){
+        for( int j = i+1; j<n; j++){
+            if(abs(mat[i][i]) < abs(mat[i][j])){
+                for(int k = 0; k<n+1;k++){
+                    //swapping matrix
+                    //			mat[i][k] = mat[i][k] + mat[j][k];
+                    //			mat[j][k] = mat[i][k] - mat[j][k];
+                    //			mat[i][k] = mat[i][k] - mat[j][k];
+                    swap(mat[i][k],mat[j][k]);
+                }
+
+            }
+        }
+    }
+    //performing gaussian elimination.
+
+    for(int i = 0; i< n-1; i++) {
+
+        for(int j = i +1; j<n;j++){
+
+            float f = mat[j][i]/mat[i][i];
+
+            for(int k = 0; k<n+1; k++){
+                mat[j][k] = mat[j][k] -f*mat[i][k];
 
-						}
-				}
-		}
-
-		//backwards substitution for unknowns.
-
-		for( int i = n-1; i >= 0; i--){
+            }
+        }
+    }
+
+    //backwards substitution for unknowns.
+
+    for( int i = n-1; i >= 0; i--){
 
-				residual[i] = mat[i][n];
-
-				for( int j = i + 1; j<n;j++){
-
-						if(i != j){
-								residual[i] = residual[i] - mat[i][j] * residual[j];
-						}
-				}
-				residual[i] = residual[i]/mat[i][i]; 
-		}
-
-
-		cout<<"A0 = "<<residual[0]<<"\n";
-		cout<<"A2 = "<<residual[1]<<"\n";
-		cout<<"A4 = "<<residual[2]<<"\n";
+        residual[i] = mat[i][n];
+
+        for( int j = i + 1; j<n;j++){
+
+            if(i != j){
+                residual[i] = residual[i] - mat[i][j] * residual[j];
+            }
+        }
+        residual[i] = residual[i]/mat[i][i]; 
+    }
+
+
+    cout<<"A0 = "<<residual[0]<<"\n";
+    cout<<"A2 = "<<residual[1]<<"\n";
+    cout<<"A4 = "<<residual[2]<<"\n";
 
 
 
 
-		if(test == 0){	
+    if(test == 0){	
 
-				printf("Please enter sigma, sigma = 0 for perfect allignment : ");
-				scanf("%lf", &Sigma);
-				if(Sigma >= 0 ){
-						sigma_token = 1; 
-				}else{
-						do{
-
-								printf("Negative numbers are not allowed!\nRe-enter Sigma : ");
-								scanf("%lf", &Sigma);
-						}while(Sigma <0);}
-
-				if(Sigma >= 0 ){
-						sigma_token = 1; printf(" --- Sigma Loaded --- \n"); 
-				}
+        printf("Please enter sigma, sigma = 0 for perfect allignment : ");
+        scanf("%lf", &Sigma);
+        if(Sigma >= 0 ){
+            sigma_token = 1; 
+        }else{
+            do{
+
+                printf("Negative numbers are not allowed!\nRe-enter Sigma : ");
+                scanf("%lf", &Sigma);
+            }while(Sigma <0);}
+
+        if(Sigma >= 0 ){
+            sigma_token = 1; printf(" --- Sigma Loaded --- \n"); 
+        }
 
-				printf("Please enter J1,J2 : ");
-				scanf("%lf,%lf",&j1,&j2);	
-				if(j1 >= 0 && j2  >= 0 ){
-						j1j2token = 1;
-				}else{
-						do{
-								printf("Negative numbers are not allowed!\nRe-enter J1,J2 : ");
-								scanf("%lf,%lf", &j1,&j2);	
-						}while(j1 < 0 || j2 < 0);}
+        printf("Please enter J1,J2 : ");
+        scanf("%lf,%lf",&j1,&j2);	
+        if(j1 >= 0 && j2  >= 0 ){
+            j1j2token = 1;
+        }else{
+            do{
+                printf("Negative numbers are not allowed!\nRe-enter J1,J2 : ");
+                scanf("%lf,%lf", &j1,&j2);	
+            }while(j1 < 0 || j2 < 0);}
 
-				if(j1 >= 0 && j2  >= 0 ){
-						j1j2token = 1; printf(" --- J1 & J2 Loaded --- \n");
-				}
+        if(j1 >= 0 && j2  >= 0 ){
+            j1j2token = 1; printf(" --- J1 & J2 Loaded --- \n");
+        }
 
-		}
+    }
 
-		// Calculate BK(j1) for perfect allignment
-		//
+    // Calculate BK(j1) for perfect allignment
+    //
 
-		double A[6];
-		double js = sqrt(2*j1 +1); 
-		double j12 = 2*j1;
-		int IS = (int)j12 % 2; 
+    double A[6];
+    double js = sqrt(2*j1 +1); 
+    double j12 = 2*j1;
+    int IS = (int)j12 % 2; 
 
-		double cg12,cg22,cg24,cg14 = 0.;
-		double I1, I2 = 0.;
+    double cg12,cg22,cg24,cg14 = 0.;
+    double I1, I2 = 0.;
 
-		double cg1, cg2 = 0.;	
+    double cg1, cg2 = 0.;	
 
-		double Bk11, Bk12 = 0.;
-		if(Sigma == 0){
-				if(IS == 1){
-						A[0] = j1;
-						A[1] = j1;
-						A[2] = 2.;
-						A[3] = .5;
-						A[4] = -.5;
-						A[5] = 0.;
+    double Bk11, Bk12 = 0.;
+    if(Sigma == 0){
+        if(IS == 1){
+            A[0] = j1;
+            A[1] = j1;
+            A[2] = 2.;
+            A[3] = .5;
+            A[4] = -.5;
+            A[5] = 0.;
 
-						cg12 = CG2(A);
+            cg12 = CG2(A);
 
-						A[3] = -.5;
-						A[4] = .5;
+            A[3] = -.5;
+            A[4] = .5;
 
-						cg22 = CG2(A);
+            cg22 = CG2(A);
 
-						A[2] = 4.;
+            A[2] = 4.;
 
-						cg24 = CG2(A);
+            cg24 = CG2(A);
 
-						A[3] = .5;
-						A[4] = -.5;
+            A[3] = .5;
+            A[4] = -.5;
 
-						cg14 = CG2(A);
+            cg14 = CG2(A);
 
 
-						printf("cg12 = %lf\n",cg12);
-						printf("cg22 = %lf\n",cg22);
-						printf("cg24 = %lf\n",cg24);
-						printf("cg24 = %lf\n",cg24);
+            printf("cg12 = %lf\n",cg12);
+            printf("cg22 = %lf\n",cg22);
+            printf("cg24 = %lf\n",cg24);
+            printf("cg24 = %lf\n",cg24);
 
-						printf("Bk11 = %lf\n",Bk11);	
-						printf("Bk12 = %lf\n",Bk12);	
+            printf("Bk11 = %lf\n",Bk11);	
+            printf("Bk12 = %lf\n",Bk12);	
 
-						Bk11 = (0.5)*js*(pow(-1,I1) * cg12 + pow(-1,I2) * cg22);
-						Bk12 = (0.5)*js*(pow(-1,I1) * cg14 + pow(-1,I2) * cg24);
-				}else if(IS == 0){
-						A[0] = j1; 
-						A[1] = j1; 
-						A[2] = 2.;
-						A[3] = 0.;
-						A[4] = 0.;
-						A[5] = 0.; 
+            Bk11 = (0.5)*js*(pow(-1,I1) * cg12 + pow(-1,I2) * cg22);
+            Bk12 = (0.5)*js*(pow(-1,I1) * cg14 + pow(-1,I2) * cg24);
+        }else if(IS == 0){
+            A[0] = j1; 
+            A[1] = j1; 
+            A[2] = 2.;
+            A[3] = 0.;
+            A[4] = 0.;
+            A[5] = 0.; 
 
-						cg1 = CG2(A);	
+            cg1 = CG2(A);	
 
-						Bk11 = pow(-1,j1) * js * cg1;
+            Bk11 = pow(-1,j1) * js * cg1;
+
+            A[2] = 4.;
 
-						A[2] = 4.;
+            cg2 = CG2(A);
+            Bk12 = pow(-1,j1) * js * cg2; 	
+
+            printf("Bk11 = %lf\n",Bk11);	
+            printf("Bk12 = %lf\n",Bk12);	
+            printf("cg1 = %lf\n",cg1);
+            printf("cg2 = %lf\n",cg2);
+        }
+    }
+    //normalize W(m1) 
+    //
+
+    j12 = 2*j1;
+    double j14 = 4*j1;
+
+    double sigsq = pow(Sigma,2);
+    double sum1 = 0.; 
+
+    double am1,amsq,x,ex = 0.;
+    //-------------------------------------------------------SOMETHING WRONG WITH CN1.
+    for(int i = 0; i <= j14; i = i + 2){
+        am1 = 0.5*(i - j12);
+        amsq = pow(am1,2);
+        x = - (amsq/(2*sigsq));
+        ex = exp(x);
+
+        sum1 = sum1 + ex;
+    }
+    double cn1 = 1./sum1;	
+
+    //	cout << "cn1"  << " " << cn1 << "\n"; 
+
+    double AL0 = 0.;
+    double A0 = j1 - j2; 
+    double A0b = abs(A0);
+    if(A0b > 0){ AL0 = A0b;}
+
+    if(A0b <= 0){AL0 = 1.;}
+
+    double AL1 = AL0 +1; 
+
+    double am11,amsq1,x1,ex1 = 0.;
+    //calculate Bk(j1) for gaussian W(m1) or non zero Sigma
+
+    //NEEDS VALUE FIX
+    if(Sigma != 0){
+        double A[6];
+        A[0] = j1; 
+        A[1] = j1; 
+        A[5] = 0.;
+        double sfact = sqrt(2*j1 +1); 
+        double cgg = 0.; 
+        double tTerm = 0.; 
+        double II = 0.;
+        for(int i = 2; i<= 4; i = i + 2){
+
+            A[2] = i; 
+  
+            for(int m = 0; m <= j14; m = m +2){
+                am11 = 0.5*(m - j12);
+                amsq1 = pow(am11,2);
+                x1 = -(amsq1/(2*sigsq));
+                II = j1 - am11; 
+                A[3] = am11; 
+                A[4] = -am11; 
+                
+                cgg = CG2(A);
+                //						cout << "cgg = " << "  " << cgg << "\n";
+                //						cout << "am11 = " << "  " << am11 << "\n";
+                //						cout << "amsq1 = " << "  " << amsq1 << "\n";
+                //						cout << "x1 = " << "  " << x1 << "\n";
+
+                ex1 = exp(x1);
+                tTerm = cn1*ex1*pow(-1,II) * sfact*cgg;
+
+                //						cout <<"tTerm "<<"  " << tTerm << "\n";
+
+                if(i == 2){ Bk11 = Bk11 + tTerm;
+                }else if(i == 4){ Bk12 = Bk12 + tTerm;} 
+            }
+        }
+        cout << "Bk11" << " " << Bk11 << "\n";
+        cout << "Bk12" << " " << Bk12 << "\n";
+
+
+    }
+
+    vector<vector<string> > content_r;
+    vector<string> row_r;
+    string line_r, word_r; 
+
+    vector<string> j1data;
+    vector<string> j2data;
+    vector<string> Rk20data; 
+    vector<string> Rk21data; 
+    vector<string> Rk22data; 
+    vector<string> Rk40data; 
+    vector<string> Rk41data; 
+    vector<string> Rk42data; 
+
+
+    fstream file1("RkTable.csv");
+    if(file1.is_open()){
+        while(getline(file1,line_r)){
+            row_r.clear();
+
+            stringstream str(line_r);
+
+            while(getline(str, word_r, ','))
+
+                row_r.push_back(word_r);
+            content_r.push_back(row_r);
+            j1data.push_back(row_r[0]);
+            j2data.push_back(row_r[1]);
+            Rk20data.push_back(row_r[2]);			
+            Rk21data.push_back(row_r[3]);			
+            Rk22data.push_back(row_r[4]);			
+            Rk40data.push_back(row_r[5]);			
+            Rk41data.push_back(row_r[6]);			
+            Rk42data.push_back(row_r[7]);			
+
+        }
+        cout<< " --- Rk(llj1j2) Table Loaded --- \n";
+    }else cout<< "Could not open the file\n";
+
+    vector<double>::iterator dit;
+    vector<double>::iterator dit2;
+
+    vector<double> j1datad   = string_to_double_vector(j1data);
+    vector<double> j2datad   = string_to_double_vector(j2data);
+    vector<double> rk20datad = string_to_double_vector(Rk20data);
+    vector<double> rk21datad = string_to_double_vector(Rk21data);
+    vector<double> rk22datad = string_to_double_vector(Rk22data);
+    vector<double> rk40datad = string_to_double_vector(Rk40data);
+    vector<double> rk41datad = string_to_double_vector(Rk41data);
+    vector<double> rk42datad = string_to_double_vector(Rk42data);
+
+    int nn = 0;
+    for(dit = j1datad.begin();dit < j1datad.end();dit++){
+
+        if(j1 == *dit){
+
+            for(dit2 = j2datad.begin();dit2 < j2datad.end();dit2++){
+
+                if(j2 == *dit2 and j1 == j1datad[distance(j2datad.begin(),dit2)]){
+                    //			cout << *dit << " " << *dit2 << " Index Equals =  " <<distance(j2datad.begin(),dit2) << "\n";
+                    nn = distance(j2datad.begin(),dit2);
+                }
+            }
+        }
+    }
+
+    //		cout << j1<< "  " << j2 << "  "<< nn << "\n";
 
-						cg2 = CG2(A);
-						Bk12 = pow(-1,j1) * js * cg2; 	
-
-						printf("Bk11 = %lf\n",Bk11);	
-						printf("Bk12 = %lf\n",Bk12);	
-						printf("cg1 = %lf\n",cg1);
-						printf("cg2 = %lf\n",cg2);
-				}
-		}
-		//normalize W(m1) 
-		//
 
-		j12 = 2*j1;
-		double j14 = 4*j1;
-
-		double sigsq = pow(Sigma,2);
-		double sum1 = 0.; 
-
-		double am1,amsq,x,ex = 0.;
-		//-------------------------------------------------------SOMETHING WRONG WITH CN1.
-		for(int i = 0; i <= j14; i = i + 2){
-				am1 = 0.5*(i - j12);
-				amsq = pow(am1,2);
-				x = - (amsq/(2*sigsq));
-				ex = exp(x);
-
-				sum1 = sum1 + ex;
-		}
-		double cn1 = 1./sum1;	
+    // rk01,rk11,rk21, rk02, rk12,rk22 //READ OUT OF Rk TABLE 
+
+    double rk01 = rk20datad[nn];
+    double rk11 = rk21datad[nn];
+    double rk21 = rk22datad[nn];
+    double rk02 = rk40datad[nn];
+    double rk12 = rk41datad[nn];
+    double rk22 = rk42datad[nn];
 
-		//	cout << "cn1"  << " " << cn1 << "\n"; 
-
-		double AL0 = 0.;
-		double A0 = j1 - j2; 
-		double A0b = abs(A0);
-		if(A0b > 0){ AL0 = A0b;}
-
-		if(A0b <= 0){AL0 = 1.;}
-
-		double AL1 = AL0 +1; 
-
-		double am11,amsq1,x1,ex1 = 0.;
-		//calculate Bk(j1) for gaussian W(m1) or non zero Sigma
-
-		//NEEDS VALUE FIX
-		if(Sigma != 0){
-				double A[6];
-				A[0] = j1; 
-				A[1] = j1; 
-
-				double sfact = sqrt(2*j1 +1); 
-				double cgg = 0.; 
-				double tTerm = 0.; 
-				double II = 0.;
-				for(int i = 2; i<= 4; i = i + 2){
 
-						A[2] = i; 
-						Bk11 = 0.;
-						Bk12 = 0.; 
-						for(int m = 0; m <= j14; m = m +2){
-								am11 = 0.5*(m - j12);
-								amsq1 = pow(am11,2);
-								x1 = -(amsq/(2*sigsq));
-								II = j1 - am11; 
-								A[3] = am11; 
-								A[4] = -am11; 
-								A[5] = 0; 
-								cgg = CG2(A);
-								//						cout << "cgg = " << "  " << cgg << "\n";
-								//						cout << "am11 = " << "  " << am11 << "\n";
-								//						cout << "amsq1 = " << "  " << amsq1 << "\n";
-								//						cout << "x1 = " << "  " << x1 << "\n";
-
-								ex1 = exp(x1);
-								tTerm = cn1*ex1*pow(-1,II) * sfact*cgg;
-
-								//						cout <<"tTerm "<<"  " << tTerm << "\n";
-
-								if(i == 2) Bk11 = Bk11 + tTerm;
-								if(i == 4) Bk12 = Bk12 + tTerm; 
-						}
-				}
-				cout << "Bk11" << " " << Bk11 << "\n";
-				cout << "Bk12" << " " << Bk12 << "\n";
-
-
-		}
-
-		vector<vector<string> > content_r;
-		vector<string> row_r;
-		string line_r, word_r; 
-
-		vector<string> j1data;
-		vector<string> j2data;
-		vector<string> Rk20data; 
-		vector<string> Rk21data; 
-		vector<string> Rk22data; 
-		vector<string> Rk40data; 
-		vector<string> Rk41data; 
-		vector<string> Rk42data; 
-
-
-		fstream file1("RkTable.csv");
-		if(file1.is_open()){
-				while(getline(file1,line_r)){
-						row_r.clear();
-
-						stringstream str(line_r);
-
-						while(getline(str, word_r, ','))
-
-								row_r.push_back(word_r);
-						content_r.push_back(row_r);
-						j1data.push_back(row_r[0]);
-						j2data.push_back(row_r[1]);
-						Rk20data.push_back(row_r[2]);			
-						Rk21data.push_back(row_r[3]);			
-						Rk22data.push_back(row_r[4]);			
-						Rk40data.push_back(row_r[5]);			
-						Rk41data.push_back(row_r[6]);			
-						Rk42data.push_back(row_r[7]);			
-
-				}
-				cout<< " --- Rk(llj1j2) Table Loaded --- \n";
-		}else cout<< "Could not open the file\n";
+    //	printf("Bk1 = %lf\n",Bk11);
+    //	printf("Bk2 = %lf\n",Bk12);
 
-		vector<double>::iterator dit;
-		vector<double>::iterator dit2;
-
-		vector<double> j1datad   = string_to_double_vector(j1data);
-		vector<double> j2datad   = string_to_double_vector(j2data);
-		vector<double> rk20datad = string_to_double_vector(Rk20data);
-		vector<double> rk21datad = string_to_double_vector(Rk21data);
-		vector<double> rk22datad = string_to_double_vector(Rk22data);
-		vector<double> rk40datad = string_to_double_vector(Rk40data);
-		vector<double> rk41datad = string_to_double_vector(Rk41data);
-		vector<double> rk42datad = string_to_double_vector(Rk42data);
-
-		int nn = 0;
-		for(dit = j1datad.begin();dit < j1datad.end();dit++){
-
-				if(j1 == *dit){
-
-						for(dit2 = j2datad.begin();dit2 < j2datad.end();dit2++){
-
-								if(j2 == *dit2 and j1 == j1datad[distance(j2datad.begin(),dit2)]){
-										//			cout << *dit << " " << *dit2 << " Index Equals =  " <<distance(j2datad.begin(),dit2) << "\n";
-										nn = distance(j2datad.begin(),dit2);
-								}
-						}
-				}
-		}
-
-		//		cout << j1<< "  " << j2 << "  "<< nn << "\n";
+    printf("rk01 = %lf\n",rk01); 
+    printf("rk11 = %lf\n",rk11); 
+    printf("rk21 = %lf\n",rk21); 
+    printf("rk02 = %lf\n",rk02); 
+    printf("rk12 = %lf\n",rk12); 
+    printf("rk22 = %lf\n",rk22); 
 
 
-		// rk01,rk11,rk21, rk02, rk12,rk22 //READ OUT OF Rk TABLE 
-
-		double rk01 = rk20datad[nn];
-		double rk11 = rk21datad[nn];
-		double rk21 = rk22datad[nn];
-		double rk02 = rk40datad[nn];
-		double rk12 = rk41datad[nn];
-		double rk22 = rk42datad[nn];
+    // here the Chi-sqs from A2 and A4 needed to be added together. 
 
 
-		//	printf("Bk1 = %lf\n",Bk11);
-		//	printf("Bk2 = %lf\n",Bk12);
+    double delta_min = -3.14159/2.;
+    double delta_max =  3.14159/2.;
+    double step = 0.001; 
 
-		printf("rk01 = %lf\n",rk01); 
-		printf("rk11 = %lf\n",rk11); 
-		printf("rk21 = %lf\n",rk21); 
-		printf("rk02 = %lf\n",rk02); 
-		printf("rk12 = %lf\n",rk12); 
-		printf("rk22 = %lf\n",rk22); 
+    double delta = 0.;
+    double atan_delta =0.;
+    double A0E = residual[0];//NEED FrOM AF FIT
+    double A2E = residual[1];
+    double A4E = residual[2]; 
+    double A2T,a2T = 0.;
+    double A4T,a4T = 0.;
+    double rd2T = 0.;
+    double rd4T = 0.;
+    double X22,X24 = 0.;
+    double X2_total = 0.;
 
+    double a2E = A2E/A0E;
+    double a4E = A4E/A0E;
 
-		// here the Chi-sqs from A2 and A4 needed to be added together. 
 
+    int points = (delta_max - delta_min) / step ; 
 
-		double delta_min = -3.14159/2.;
-		double delta_max =  3.14159/2.;
-		double step = 0.001; 
+    vector<double> chisqr;
+    vector<double> tdelta;
 
-		double delta = 0.;
-		double atan_delta =0.;
-		double A0E = residual[0];//NEED FrOM AF FIT
-		double A2E = residual[1];
-		double A4E = residual[2]; 
-		double A2T,a2T = 0.;
-		double A4T,a4T = 0.;
-		double rd2T = 0.;
-		double rd4T = 0.;
-		double X22,X24 = 0.;
-		double X2_total = 0.;
+    for(int i = 0; i< points; i++){
+        atan_delta = i*step + delta_min;
+        delta = tan(atan_delta);
+        rd2T = (rk01 + 2*delta*rk11 + pow(delta,2)*rk21)/(1+pow(delta,2));
+        A2T = QD2*Bk11*rd2T;
+        a2T = A2T/A0E;
+        //using the idea that A0E is our normalizer. 
 
-		double a2E = A2E/A0E;
-		double a4E = A4E/A0E;
+        rd4T = (rk02 + 2*delta*rk12 + pow(delta,2)*rk22)/(1+pow(delta,2));
+        A4T = QD4*Bk12*rd4T;
+        a4T = A4T/A0E;
 
+        X22 = pow(abs(a2E -a2T),2)/(abs(a2T));
+        X24 = pow(abs(a4E -a4T),2)/(abs(a4T));
+        X2_total = (X22 + X24)/2;
 
-		int points = (delta_max - delta_min) / step ; 
+        chisqr.push_back(-log(X2_total));
+        tdelta.push_back(atan_delta);
 
-		vector<double> chisqr;
-		vector<double> tdelta;
+    }
 
-		for(int i = 0; i< points; i++){
-				atan_delta = i*step + delta_min;
-				delta = tan(atan_delta);
-				rd2T = (rk01 + 2*delta*rk11 + pow(delta,2)*rk21)/(1+pow(delta,2));
-				A2T = QD2*Bk11*rd2T;
-				a2T = A2T/A0E;
-				//using the idea that A0E is our normalizer. 
+    vector<double> Theta;
+    vector<double> AD_I;
+    double ad_start = 0.;	
+    int adpoints = (3.1415)/(step);
+    double theta,Iad = 0.;
+    double aaa,aab,aac = 0.;
+    for(int i = 0; i < adpoints; i++){
 
-				rd4T = (rk02 + 2*delta*rk12 + pow(delta,2)*rk22)/(1+pow(delta,2));
-				A4T = QD4*Bk12*rd4T;
-				a4T = A4T/A0E;
+        theta = i*step + ad_start;
+        aaa = 1.;
+        aab  = (A2E/A0E)*(1.5 * pow(cos(theta),2) - .5);
+        aac  = (A4E/A0E)*(35./8. * pow(cos(theta),4) - 30./8. * pow(cos(theta),2)  +  3./8. );
 
-				X22 = pow(abs(a2E -a2T),2)/(abs(a2T));
-				X24 = pow(abs(a4E -a4T),2)/(abs(a4T));
-				X2_total = (X22 + X24)/2;
+        Iad = aaa + aab + aac;
 
-				chisqr.push_back(-log(X2_total));
-				tdelta.push_back(atan_delta);
+        AD_I.push_back(Iad);
+        Theta.push_back(theta);
+        //these arrays arent used, but they are the function (if one just wanted to the legendre fit function that goes over the data. 
+    }
 
-		}
 
-		vector<double> Theta;
-		vector<double> AD_I;
-		double ad_start = 0.;	
-		int adpoints = (3.1415)/(step);
-		double theta,Iad = 0.;
-		double aaa,aab,aac = 0.;
-		for(int i = 0; i < adpoints; i++){
+    vector<double> dydatas;
 
-				theta = i*step + ad_start;
-				aaa = 1.;
-				aab  = (A2E/A0E)*(1.5 * pow(cos(theta),2) - .5);
-				aac  = (A4E/A0E)*(35./8. * pow(cos(theta),4) - 30./8. * pow(cos(theta),2)  +  3./8. );
+    for(int i = 0; i < dydata.size(); i++){
+        double bbb = dydata[i]; 
 
-				Iad = aaa + aab + aac;
+        dydatas.push_back(bbb/A0E);
+        printf("Normalized y-intensity %i = %lf\n",i+1,dydatas[i]);
 
-				AD_I.push_back(Iad);
-				Theta.push_back(theta);
-				//these arrays arent used, but they are the function (if one just wanted to the legendre fit function that goes over the data. 
-		}
+    }
+    //This is to temporarily fix the issue with the AD distribution, it takes off the last point in the arrays when plotting or doesnt see them. NEED FIX
+    dydatas.push_back(A0E);
+    dangler.push_back(3.1415/2.);
 
+    int param;
 
-		vector<double> dydatas;
+    param = param_run(det_param_token, gamma_energy_token, ang_file_token, sigma_token, j1j2token);  	
+    /*	
+        cout << det_param_token << "\n"; 
+        cout << gamma_energy_token <<  "\n";  
+        cout << ang_file_token <<  "\n"; 
+        cout << sigma_token <<  "\n";  
+        cout << j1j2token <<  "\n"; 
 
-		for(int i = 0; i < dydata.size(); i++){
-				double bbb = dydata[i]; 
+        cout << param << "\n";
+        */
 
-				dydatas.push_back(bbb/A0E);
-				printf("Normalized y-intensity %i = %lf\n",i+1,dydatas[i]);
+    //Initialize Graphics Here. 
 
-		}
-		//This is to temporarily fix the issue with the AD distribution, it takes off the last point in the arrays when plotting or doesnt see them. NEED FIX
-		dydatas.push_back(A0E);
-		dangler.push_back(3.1415/2.);
 
-		int param;
+    HistoGUI gui;
 
-		param = param_run(det_param_token, gamma_energy_token, ang_file_token, sigma_token, j1j2token);  	
-		/*	
-			cout << det_param_token << "\n"; 
-			cout << gamma_energy_token <<  "\n";  
-			cout << ang_file_token <<  "\n"; 
-			cout << sigma_token <<  "\n";  
-			cout << j1j2token <<  "\n"; 
+    HistoGUIad gui_ad; 
 
-			cout << param << "\n";
-		 */
+// FIX menu input complexity. 
+//
+    int optnum = -1;
+    menu();
 
-		//Initialize Graphics Here. 
+    printf("Input option : ");
+    scanf("%d", &optnum);
+    if(optnum < 0){
+        optnum = -1;
+        printf("Negative numbers are not allowed!\nInput option : ");
+        scanf("%d", &optnum);
+    }
 
+    if(optnum == 0){
+        Readme();
+        optnum = -1;
 
-		HistoGUI gui;
+        printf("Input option : ");
+        scanf("%d", &optnum);
+    } 
 
-		HistoGUIad gui_ad; 
+    //if all inputs are valid, plot distribution. 
+    if(param == 1 && optnum == 1){
+        optnum = -1;
 
 
-		int optnum = -1;
-		menu();
+        //plotting values here
+        //x will be arctan of mixing ratio. 
+        //y will be log(X^2); 
 
-		printf("Input option : ");
-		scanf("%d", &optnum);
-		if(optnum < 0){
-				optnum = -1;
-				printf("Negative numbers are not allowed!\nInput option : ");
-				scanf("%d", &optnum);
-		}
+        gui.SetData(tdelta,chisqr);
+        //gui.SetData(Theta,AD_I);
+        gui.Init();
+        gui.Loop();
+        gui.Close();
 
-		if(optnum == 0){
-				Readme();
-				optnum = -1;
+        //		printf("Input option : ");
+        //		scanf("%d", &optnum);
+        //		if(optnum < 0){		
+        //			printf("Input option : ");
+        //			scanf("%d", &optnum);
+        //		}
+    }
 
-				printf("Input option : ");
-				scanf("%d", &optnum);
-		} 
+    //Plot Angular distribution fit, with points and error bars. 
+    if(param == 1 and optnum == 2){
+        /*
+           gui.SetData(Theta,AD_I);
+           gui.Init();
+           gui.Loop();
+           gui.Close();
+           */
 
-		//if all inputs are valid, plot distribution. 
-		if(param == 1 && optnum == 1){
-				optnum = -1;
+        gui_ad.SetData(dangler,dydatas);
+        gui_ad.SetErrors(deydata);
+        gui_ad.SetFit(residual[0],residual[1],residual[2]);
 
+        gui_ad.Init();
+        gui_ad.Loop();
+        gui_ad.Close();
 
-				//plotting values here
-				//x will be arctan of mixing ratio. 
-				//y will be log(X^2); 
 
-				gui.SetData(tdelta,chisqr);
-				//gui.SetData(Theta,AD_I);
-				gui.Init();
-				gui.Loop();
-				gui.Close();
 
-				//		printf("Input option : ");
-				//		scanf("%d", &optnum);
-				//		if(optnum < 0){		
-				//			printf("Input option : ");
-				//			scanf("%d", &optnum);
-				//		}
-		}
+    }
+    //exit
 
-		//Plot Angular distribution fit, with points and error bars. 
-		if(param == 1 and optnum == 2){
-				/*
-				   gui.SetData(Theta,AD_I);
-				   gui.Init();
-				   gui.Loop();
-				   gui.Close();
-				 */
+    printf("Input option : ");
+    scanf("%d", &optnum);
+    if(optnum < 0){
+        optnum = -1;
+        printf("Negative numbers are not allowed!\nInput option : ");
+        scanf("%d", &optnum);
+    }
 
-				gui_ad.SetData(dangler,dydatas);
-				gui_ad.SetErrors(deydata);
-				gui_ad.SetFit(residual[0],residual[1],residual[2]);
+    if(param == 1 and optnum == 3){
 
-				gui_ad.Init();
-				gui_ad.Loop();
-				gui_ad.Close();
+        return 0;
+    }
 
-
-
-		}
-		//exit
-		if(param == 1 and optnum == 3){
-
-				return 0;
-		}
-
-
-		return 1;
+    return 1;
 
 
 }
diff --git a/GUI.h b/Graphical_Testing/GUI_fit_test.h
similarity index 100%
rename from GUI.h
rename to Graphical_Testing/GUI_fit_test.h
diff --git a/LiveLook_2d.cxx b/Graphical_Testing/LiveLook_2d.cxx
similarity index 100%
rename from LiveLook_2d.cxx
rename to Graphical_Testing/LiveLook_2d.cxx
diff --git a/Racah.h b/Racah.h
index d354b95..f825db9 100644
--- a/Racah.h
+++ b/Racah.h
@@ -5,7 +5,7 @@
 
 using namespace std;
 
-//this program doesnt work and it is not used, I spent so much time transcribing it that maybe, one day some poor soul will come and try to fix it to match rose and brinks. 
+//this program doesnt work and it is not used, I spent so much time transcribing it that maybe, one day some poor soul will come and try to fix it to match the table in rose and brinks. 
 
 double RACAH(double A[6]){
 	
diff --git a/Test_Scripts/bk_test.cxx b/Test_Scripts/bk_test.cxx
new file mode 100644
index 0000000..5bb2205
--- /dev/null
+++ b/Test_Scripts/bk_test.cxx
@@ -0,0 +1,117 @@
+#include "global.h"
+using namespace std;
+int main(int argc, char** argv){
+
+    double j1 = 2.;
+    double j2 = 1.;
+    double Sigma = 0.5;
+
+    double j12 = 2*j1;
+    double j14 = 4*j1;
+
+    double sigsq = pow(Sigma,2);
+    double sum1 = 0.; 
+
+    double am1,amsq,x,ex = 0.;
+    //-------------------SOMETHING WRONG WITH CN1.
+    for(int i = 0; i <= j14; i = i + 2){
+        am1 = 0.5*(i - j12);
+        amsq = pow(am1,2);
+        x = - (amsq/(2*sigsq));
+        ex = exp(x);
+
+        sum1 = sum1 + ex;
+ //       cout << "sum 1 = " << sum1 << "\n";
+    }
+    double cn1 = 1./sum1;	
+
+    cout << "cn1"  << " " << cn1 << "\n"; 
+    
+    double am11,amsq1,x1,ex1 = 0.;
+    //calculate Bk(j1) for gaussian W(m1) or non zero Sigma
+
+    double Bk11,Bk12 = 0.;
+
+    //NEEDS VALUE FIX
+    if(Sigma != 0){
+        double A[6];
+        A[0] = j1; 
+        A[1] = j1; 
+        A[5] = 0;
+
+        double sfact = sqrt(2*j1 +1); 
+        double cgg = 0.; 
+        double tTerm = 0.; 
+        double II = 0.;
+        for(int i = 2; i<= 4; i = i + 2){
+
+            A[2] = i; 
+  
+            for(int m = 0; m <= j14; m = m +2){
+                am11 = 0.5*(m - j12);
+                amsq1 = pow(am11,2);
+
+                x1 = -(amsq1/(2*sigsq));
+                II = j1 - am11; 
+                A[3] = am11; 
+                A[4] = -am11; 
+/*           
+                cout << "---------\n";
+                cout << "i = " << i << "\n";
+                cout << "m = " << m << "\n";
+                
+                cout << "AM1 = " << am11 << "\n";       
+                cout << "AMSQ = " << amsq1 << "\n";       
+                cout << "X = " << x1 << "\n";       
+                cout << "I = " << II << "\n";       
+                
+       
+         
+                cout << "---------\n";
+*/         
+         
+                cgg = 0.5;  
+
+                ex1 = exp(x1);
+                tTerm = cn1*ex1*pow(-1,II) * sfact*cgg;
+/*        
+                cout << "---------\n";
+                cout << "i = " << i << "\n";
+                cout << "m = " << m << "\n";
+                cout << "ex1 = " << ex1 << "\n";
+                cout << "tTerm = " << tTerm << "\n";
+*/
+                cout << "---------\n";
+                if(i == 2){ 
+                
+                Bk11 = Bk11 + tTerm;
+                cout << "---------\n";
+                cout << "i = " << i << "\n";
+                cout << "m = " << m << "\n";
+                cout << "tTerm = " << tTerm << "\n";
+                cout << "BK1 =" << Bk11 << "\n";
+                cout << "---------\n";
+                }
+                else if(i == 4){ 
+                Bk12 = Bk12 + tTerm; 
+                
+                cout << "---------\n";
+                cout << "i = " << i << "\n";
+                cout << "m = " << m << "\n";
+                cout << "tTerm = " << tTerm << "\n";
+                cout << "BK2 =" << Bk12 << "\n";
+                cout << "---------\n";
+
+
+                }
+            }
+        }
+        cout << "Bk11" << " " << Bk11 << "\n";
+        cout << "Bk12" << " " << Bk12 << "\n";
+
+
+    }
+
+    return 0; 
+
+}
diff --git a/racah_table_reader_test.cxx b/Test_Scripts/racah_table_reader_test.cxx
similarity index 100%
rename from racah_table_reader_test.cxx
rename to Test_Scripts/racah_table_reader_test.cxx