-// This function is taken from Numerical Recipes in C++ by Press et al., 2nd edition, pg. 479
-int MothurMetastats::spline(vector<double>& x, vector<double>& y, int yp1, int ypn, vector<double>& y2) {
- try {
-
- cout << "lambdas" << endl;
- for (int l = 0; l < x.size(); l++){ cout << x[l] << '\t'; }
- cout << endl << "pi0_hat" << endl;
- for (int l = 0; l < y.size(); l++){ cout << y[l] << '\t'; }
- cout << endl;
-
- double p, qn, sig, un;
-
- int n = y2.size();
- vector<double> u(n-1, 0.0);
-
- if (yp1 > 0.99e30) { y2[0] = u[0] = 0.0; }
- else {
- y2[0] = -0.5;
- u[0] = (3.0 / (x[1] - x[0])) * ((y[1] - y[0]) / (x[1]-x[0]) - yp1);
- }
-
- for (int i = 1; i < n-1; i++) {
- sig = (x[i]-x[i-1])/(x[i+1]-x[i-1]);
- p = sig * y2[i-1] + 2.0;
- y2[i] = (sig - 1.0)/p;
- u[i] = (y[i+1]-y[i]) / (x[i+1]-x[i]) - (y[i] - y[i-1]) / (x[i]-x[i-1]);
- u[i] = (6.0*u[i]/(x[i+1]-x[i-1])-sig*u[i-1])/p;
- }
-
- if (ypn > 0.99e30) { qn=un=0.0; }
- else {
- qn=0.5;
- un=(3.0/(x[n-1]-x[n-2]))*(ypn-(y[n-1]-y[n-2])/(x[n-1]-x[n-2]));
- }
-
- y2[n-1]=(un-qn*u[n-2])/(qn*y2[n-2]+1.0);
-
- for (int k=n-2; k>=0; k--) {
- y2[k]=y2[k]*y2[k+1]+u[k];
- }
-
- return 0;
-
- }catch(exception& e) {
- m->errorOut(e, "MothurMetastats", "spline");
- exit(1);
- }
-}
-/***********************************************************/
-// This function is taken from Numerical Recipes in C++ by Press et al., 2nd edition, pg. 479
-int MothurMetastats::splint(vector<double>& xa, vector<double>& ya, double& x, double& y, vector<double>& y2a) {
- try {
- int k;
- double h,b,a;
-
- int n = xa.size();
-
- int klo=0;
- int khi=n-1;
- while (khi-klo > 1) {
-
- if (m->control_pressed) { break; }
-
- k = (khi+klo) >> 1;
- if (xa[k] > x) { khi=k; }
- else { klo=k; }
- }
-
- h=xa[khi]-xa[klo];
- if (h == 0.0) { m->mothurOut("[ERROR]: Bad xa input to splint routine."); m->mothurOutEndLine(); m->control_pressed = true; return 0; }
- a=(xa[khi]-x)/h;
- b=(x - xa[klo])/h;
- y=a*ya[klo]+b*ya[khi]+((a*a*a-a)*y2a[klo]+(b*b*b-b)*y2a[khi])*(h*h)/6.0;
-
- return 0;
-
- }catch(exception& e) {
- m->errorOut(e, "MothurMetastats", "splint");
- exit(1);
- }
-}
-/***********************************************************/
-vector<double> MothurMetastats::sknot1(vector<double> x) {