/******************************************************************************/
/*									      */
/*  COMPUTE UPPER BOUNDS OF MINIMA					      */
/*  08/25/10 (dkc)							      */
/*									      */
/*  This C program computes upper bounds of minima where l is prime.  This    */
/*  program is for use on the C64.					      */
/*									      */
/******************************************************************************/
#include <math.h>
void mul3shft(unsigned int *c, unsigned int b, unsigned int n);
void setn(unsigned int *a, unsigned int b, unsigned int n);

far unsigned int XG[198140];	// global copy of X
far double max0[400],min0[400]; // maximum and minimum upper bounds of minima
far unsigned int maxi[400];	// maximum upper bounds of minima (integer)
far unsigned int output[400];	// check for overflow of X array
far unsigned int hist0[400];	// histogram of x values
far unsigned int error[6];	// error array
far double outx[2];		// maximum and minimum upper bounds of minima
far unsigned int outp[2];	// maximum and minimum l values
far unsigned int table[667]={	// prime look-up table
 7, 11, 13, 17, 19, 23, 29, 31,
 37, 41, 43, 47, 53, 59, 61, 67, 71, 73,
 79, 83, 89, 97, 101, 103, 107, 109, 113, 127,
 131, 137, 139, 149, 151, 157, 163, 167, 173, 179,
 181, 191, 193, 197, 199, 211, 223, 227, 229, 233,
 239, 241, 251, 257, 263, 269, 271, 277, 281, 283,
 293, 307, 311, 313, 317, 331, 337, 347, 349, 353,
 359, 367, 373, 379, 383, 389, 397, 401, 409, 419,
 421, 431, 433, 439, 443, 449, 457, 461, 463, 467,
 479, 487, 491, 499, 503, 509, 521, 523, 541, 547,
 557, 563, 569, 571, 577, 587, 593, 599, 601, 607,
 613, 617, 619, 631, 641, 643, 647, 653, 659, 661,
 673, 677, 683, 691, 701, 709, 719, 727, 733, 739,
 743, 751, 757, 761, 769, 773, 787, 797, 809, 811,
 821, 823, 827, 829, 839, 853, 857, 859, 863, 877,
 881, 883, 887, 907, 911, 919, 929, 937, 941, 947,
 953, 967, 971, 977, 983, 991, 997, 1009, 1013, 1019,
 1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063, 1069, 1087,
 1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129, 1151, 1153,
 1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217, 1223, 1229,
 1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289, 1291, 1297,
 1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367, 1373, 1381,
 1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447, 1451, 1453,
 1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499, 1511, 1523,
 1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579, 1583, 1597,
 1601, 1607, 1609, 1613, 1619, 1621, 1627, 1637, 1657, 1663,
 1667, 1669, 1693, 1697, 1699, 1709, 1721, 1723, 1733, 1741,
 1747, 1753, 1759, 1777, 1783, 1787, 1789, 1801, 1811, 1823,
 1831, 1847, 1861, 1867, 1871, 1873, 1877, 1879, 1889, 1901,
 1907, 1913, 1931, 1933, 1949, 1951, 1973, 1979, 1987, 1993,
 1997, 1999, 2003, 2011, 2017, 2027, 2029, 2039, 2053, 2063,
 2069, 2081, 2083, 2087, 2089, 2099, 2111, 2113, 2129, 2131,
 2137, 2141, 2143, 2153, 2161, 2179, 2203, 2207, 2213, 2221,
 2237, 2239, 2243, 2251, 2267, 2269, 2273, 2281, 2287, 2293,
 2297, 2309, 2311, 2333, 2339, 2341, 2347, 2351, 2357, 2371,
 2377, 2381, 2383, 2389, 2393, 2399, 2411, 2417, 2423, 2437,
 2441, 2447, 2459, 2467, 2473, 2477, 2503, 2521, 2531, 2539,
 2543, 2549, 2551, 2557, 2579, 2591, 2593, 2609, 2617, 2621,
 2633, 2647, 2657, 2659, 2663, 2671, 2677, 2683, 2687, 2689,
 2693, 2699, 2707, 2711, 2713, 2719, 2729, 2731, 2741, 2749,
 2753, 2767, 2777, 2789, 2791, 2797, 2801, 2803, 2819, 2833,
 2837, 2843, 2851, 2857, 2861, 2879, 2887, 2897, 2903, 2909,
 2917, 2927, 2939, 2953, 2957, 2963, 2969, 2971, 2999, 3001,
 3011, 3019, 3023, 3037, 3041, 3049, 3061, 3067, 3079, 3083,
 3089, 3109, 3119, 3121, 3137, 3163, 3167, 3169, 3181, 3187,
 3191, 3203, 3209, 3217, 3221, 3229, 3251, 3253, 3257, 3259,
 3271, 3299, 3301, 3307, 3313, 3319, 3323, 3329, 3331, 3343,
 3347, 3359, 3361, 3371, 3373, 3389, 3391, 3407, 3413, 3433,
 3449, 3457, 3461, 3463, 3467, 3469, 3491, 3499, 3511, 3517,
 3527, 3529, 3533, 3539, 3541, 3547, 3557, 3559, 3571, 3581,
 3583, 3593, 3607, 3613, 3617, 3623, 3631, 3637, 3643, 3659,
 3671, 3673, 3677, 3691, 3697, 3701, 3709, 3719, 3727, 3733,
 3739, 3761, 3767, 3769, 3779, 3793, 3797, 3803, 3821, 3823,
 3833, 3847, 3851, 3853, 3863, 3877, 3881, 3889, 3907, 3911,
 3917, 3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989, 4001,
 4003, 4007, 4013, 4019, 4021, 4027, 4049, 4051, 4057, 4073,
 4079, 4091, 4093, 4099, 4111, 4127, 4129, 4133, 4139, 4153,
 4157, 4159, 4177, 4201, 4211, 4217, 4219, 4229, 4231, 4241,
 4243, 4253, 4259, 4261, 4271, 4273, 4283, 4289, 4297, 4327,
 4337, 4339, 4349, 4357, 4363, 4373, 4391, 4397, 4409, 4421,
 4423, 4441, 4447, 4451, 4457, 4463, 4481, 4483, 4493, 4507,
 4513, 4517, 4519, 4523, 4547, 4549, 4561, 4567, 4583, 4591,
 4597, 4603, 4621, 4637, 4639, 4643, 4649, 4651, 4657, 4663,
 4673, 4679, 4691, 4703, 4721, 4723, 4729, 4733, 4751, 4759,
 4783, 4787, 4789, 4793, 4799, 4801, 4813, 4817, 4831, 4861,
 4871, 4877, 4889, 4903, 4909, 4919, 4931, 4933, 4937, 4943,
 4951, 4957, 4967, 4969, 4973, 4987, 4993, 4999, 5003};

int main () {
unsigned int maxn=4000000;     // maximum n value
unsigned int i,j,k,a,b,p,maxp,minp;
double x,ln2,ln3,maxx,minx;
double x0,x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11,x12,x13,x14,x15,x16,x17;
unsigned int iter0,iter1,iter2,iter3,iter4,iter5,iter6,iter7,iter8,iter9,iter10;
unsigned int iter11,iter12,iter13,iter14,iter15,iter16,iter17;
unsigned int save0,save1,save2,save3,save4,save5,save6,save7,save8,save9,save10;
unsigned int save11,save12,save13,save14,save15,save16,save17;
unsigned int last0,last1,last2,last3,last4,last5,last6,last7,last8,last9,last10;
unsigned int last11,last12,last13,last14,last15,last16,last17;
unsigned int delta,count,limit,n;
unsigned int X[198140];
if (maxn<=165392)		// set number of words
   n=8192;
else {
   if (maxn<=330784)
      n=16384;
   else {
      if (maxn<=661576)
	 n=32768;
      else {
	 if (maxn<=1323152)
	    n=65536;
	 else {
	    if (maxn<=2646304)
	       n=131072;
	    else {
	       if (maxn<=4000000)
		  n=198140;
	       else{
		  error[0]=1;	// not enough internal memory
		  n=198140;
		  }
	       }
	    }
	 }
      }
   }
for (i=0; i<2; i++) {		// clear maximum/minimum array
   outx[i]=0.0;
   outp[i]=0;
   }
maxx=-655360.0; 		// set maximum
minx=65536.0;			// set minimum
for (i=0; i<400; i++) { 	// clear histogram
   hist0[i]=0;
   output[i]=0;
   max0[i]=0.0;
   min0[i]=65536.0*65536.0*65536.0*65536.0;
   maxi[i]=0;
   }
count=0;			// clear counter
error[0]=0;			// clear error array
error[1]=0;
error[2]=0;
error[3]=0;
error[4]=0;
error[5]=0;
iter0=0;			// period count
iter1=0;
iter2=0;
iter3=0;
iter4=0;
iter5=0;
iter6=0;
iter7=0;
iter8=0;
iter9=0;
iter10=0;
iter11=0;
iter12=0;
iter13=0;
iter14=0;
iter15=0;
iter16=0;
iter17=0;
save0=0;			// prime from last period
save1=0;
save2=0;
save3=0;
save4=0;
save5=0;
save6=0;
save7=0;
save8=0;
save9=0;
save10=0;
save11=0;
save12=0;
save13=0;
save14=0;
save15=0;
save16=0;
save17=0;
last0=0;			// last prime
last1=0;
last2=0;
last3=0;
last4=0;
last5=0;
last6=0;
last7=0;
last8=0;
last9=0;
last10=0;
last11=0;
last12=0;
last13=0;
last14=0;
last15=0;
last16=0;
last17=0;
x0=0.0; 			// last minimum
x1=0.0;
x2=0.0;
x3=0.0;
x4=0.0;
x5=0.0;
x6=0.0;
x7=0.0;
x8=0.0;
x9=0.0;
x10=0.0;
x11=0.0;
x12=0.0;
x13=0.0;
x14=0.0;
x15=0.0;
x16=0.0;
x17=0.0;
ln2=log(2);			// log(2)
ln3=log(3);			// log(3)
a=0;				// clear a
b=0;				// clear b
setn(X, 0, n);
X[0]=0x10000000;		// x=1/2
for (i=0; i<maxn; i++) {	// maximum n value
   if (((X[0]+0xe0000000)&0x80000000)==0) {  // x>1
      mul3shft(X, 2, n);	// x=(3/4)x
      a=a+1;
      }
   else {			// x<1
      mul3shft(X, 1, n);	// x=(3/2)x
      b=b+1;
      }
   if (i>=(maxn-400)) { 	// check for overflow
      output[count]=X[n-10];
      count=count+1;
      }
   error[2]=a;
   error[3]=b;
   j=2*a+b+1;			   // l
   if (j==1)			   // check if l=1
      goto bskip;
   if ((j&1)==0)		   // check if 2 divides l
      goto bskip;
   if (j==(j/3)*3)		   // check if 3 divides l
      goto bskip;
   if (j==(j/5)*5)		   // check if 5 divides l
      goto bskip;
   if (j>table[666])		   // check if l is greater than last prime
      goto yskip;
   for (k=0; k<667; k++) {	   // check if l is in the prime LUT
      if (j==table[k])
	 goto askip;
      if (table[k]>j)
	 goto bskip;
      }
   goto bskip;
yskip:
   limit=(unsigned int)(sqrt((double)j)+20.0);
   for (k=0; k<667; k++) {
      p=table[k];		   // load prime
      if (p>limit)		   // check limit
	 goto askip;
      if (j==(j/p)*p)		   // check if l is divisible by prime
	 goto bskip;
      }
askip:
   error[4]=j;			   // save l
   x=(double)(a+1)/((double)(2*a+b+1)*ln2-(double)(a+b)*ln3);  // upper bound of minimum
   if (x>maxx) {		   // compare to maximum
      maxx=x;
      maxp=j;
      }
   if (x<minx) {		   // compare to minimum
      minx=x;
      minp=j;
      }
   if ((j+1)==((j+1)/54)*54) {	   // check if 54 divides l+1
      if ((x0>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last0-save0)/54;   // difference in l values
	 save0=last0;		   // save old l value
	 if (iter0!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x0)
	       max0[delta]=x0;
	    if (min0[delta]>x0)
	       min0[delta]=x0;
	    }
	 iter0=iter0+1; 	   // increment period count
	 }
      last0=j;			   // save l value
      x0=x;			   // last x value
      }
   if ((j+5)==((j+5)/54)*54) {	   // check if 54 divides l+5
      if ((x1>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last1-save1)/54;   // difference in l values
	 save1=last1;		   // save old l value
	 if (iter1!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x1)
	       max0[delta]=x1;
	    if (min0[delta]>x1)
	       min0[delta]=x1;
	    }
	 iter1=iter1+1; 	   // increment period count
	 }
      last1=j;			   // save l value
      x1=x;			   // last x value
      }
   if ((j+7)==((j+7)/54)*54) {	   // check if 54 divides l+7
      if ((x2>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last2-save2)/54;   // difference in l values
	 save2=last2;		   // save old l value
	 if (iter2!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x2)
	       max0[delta]=x2;
	    if (min0[delta]>x2)
	       min0[delta]=x2;
	    }
	 iter2=iter2+1; 	   // increment period count
	 }
      last2=j;			   // save l value
      x2=x;			   // last x value
      }
   if ((j+11)==((j+11)/54)*54) {   // check if 54 divides l+11
      if ((x3>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last3-save3)/54;   // difference in l values
	 save3=last3;		   // save old l value
	 if (iter3!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x3)
	       max0[delta]=x3;
	    if (min0[delta]>x3)
	       min0[delta]=x3;
	    }
	 iter3=iter3+1; 	   // increment period count
	 }
      last3=j;			   // save l value
      x3=x;			   // last x value
      }
   if ((j+13)==((j+13)/54)*54) {   // check if 54 divides l+13
      if ((x4>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last4-save4)/54;   // difference in l values
	 save4=last4;		   // save old l value
	 if (iter4!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x4)
	       max0[delta]=x4;
	    if (min0[delta]>x4)
	       min0[delta]=x4;
	    }
	 iter4=iter4+1; 	   // increment period count
	 }
      last4=j;			   // save l value
      x4=x;			   // last x value
      }
   if ((j+17)==((j+17)/54)*54) {   // check if 54 divides l+17
      if ((x5>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last5-save5)/54;   // difference in l values
	 save5=last5;		   // save old l value
	 if (iter5!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x5)
	       max0[delta]=x5;
	    if (min0[delta]>x5)
	       min0[delta]=x5;
	    }
	 iter5=iter5+1; 	   // increment period count
	 }
      last5=j;			   // save l value
      x5=x;			   // last x value
      }
   if ((j+19)==((j+19)/54)*54) {   // check if 54 divides l+19
      if ((x6>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last6-save6)/54;   // difference in l values
	 save6=last6;		   // save old l value
	 if (iter6!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x6)
	       max0[delta]=x6;
	    if (min0[delta]>x6)
	       min0[delta]=x6;
	    }
	 iter6=iter6+1; 	   // increment period count
	 }
      last6=j;			   // save l value
      x6=x;			   // last x value
      }
   if ((j+23)==((j+23)/54)*54) {   // check if 54 divides l+23
      if ((x7>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last7-save7)/54;   // difference in l values
	 save7=last7;		   // save old l value
	 if (iter7!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x7)
	       max0[delta]=x7;
	    if (min0[delta]>x7)
	       min0[delta]=x7;
	    }
	 iter7=iter7+1; 	   // increment period count
	 }
      last7=j;			   // save l value
      x7=x;			   // last x value
      }
   if ((j+25)==((j+25)/54)*54) {   // check if 54 divides l+25
      if ((x8>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last8-save8)/54;   // difference in l values
	 save8=last8;		   // save old l value
	 if (iter8!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x8)
	       max0[delta]=x8;
	    if (min0[delta]>x8)
	       min0[delta]=x8;
	    }
	 iter8=iter8+1; 	   // increment period count
	 }
      last8=j;			   // save l value
      x8=x;			   // last x value
      }
   if ((j+29)==((j+29)/54)*54) {   // check if 54 divides n+29
      if ((x9>0.0)&&(x<0.0)) {	   // check for new period
	 delta=(last9-save9)/54;   // difference in l values
	 save9=last9;		   // save old l value
	 if (iter9!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x9)
	       max0[delta]=x9;
	    if (min0[delta]>x9)
	       min0[delta]=x9;
	    }
	 iter9=iter9+1; 	   // increment period count
	 }
      last9=j;			   // save l value
      x9=x;			   // last x value
      }
   if ((j+31)==((j+31)/54)*54) {   // check if 54 divides n+31
      if ((x10>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last10-save10)/54; // difference in l values
	 save10=last10; 	   // save old l value
	 if (iter10!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x10)
	       max0[delta]=x10;
	    if (min0[delta]>x10)
	       min0[delta]=x10;
	    }
	 iter10=iter10+1;	   // increment period count
	 }
      last10=j; 		   // save l value
      x10=x;			   // last x value
      }
   if ((j+35)==((j+35)/54)*54) {   // check if 54 divides l+35
      if ((x11>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last11-save11)/54; // difference in l values
	 save11=last11; 	   // save old l value
	 if (iter11!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x11)
	       max0[delta]=x11;
	    if (min0[delta]>x11)
	       min0[delta]=x11;
	    }
	 iter11=iter11+1;	   // increment period count
	 }
      last11=j; 		   // save l value
      x11=x;			   // last x value
      }
   if ((j+37)==((j+37)/54)*54) {   // check if 54 divides l+37
      if ((x12>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last12-save12)/54; // difference in l values
	 save12=last12; 	   // save old l value
	 if (iter12!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x12)
	       max0[delta]=x12;
	    if (min0[delta]>x12)
	       min0[delta]=x12;
	    }
	 iter12=iter12+1;	   // increment period count
	 }
      last12=j; 		   // save l value
      x12=x;			   // last x value
      }
   if ((j+41)==((j+41)/54)*54) {   // check if 54 divides l+41
      if ((x13>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last13-save13)/54; // difference in l values
	 save13=last13; 	   // save old l value
	 if (iter13!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x13)
	       max0[delta]=x13;
	    if (min0[delta]>x13)
	       min0[delta]=x13;
	    }
	 iter13=iter13+1;	   // increment period count
	 }
      last13=j; 		   // save l value
      x13=x;			   // last x value
      }
   if ((j+43)==((j+43)/54)*54) {   // check if 54 divides l+43
      if ((x14>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last14-save14)/54; // difference in l values
	 save14=last14; 	   // save old l value
	 if (iter14!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x14)
	       max0[delta]=x14;
	    if (min0[delta]>x14)
	       min0[delta]=x14;
	    }
	 iter14=iter14+1;	   // increment period count
	 }
      last14=j; 		   // save l value
      x14=x;			   // last x value
      }
   if ((j+47)==((j+47)/54)*54) {   // check if 54 divides l+47
      if ((x15>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last15-save15)/54; // difference in l values
	 save15=last15; 	   // save old l value
	 if (iter15!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x15)
	       max0[delta]=x15;
	    if (min0[delta]>x15)
	       min0[delta]=x15;
	    }
	 iter15=iter15+1;	   // increment period count
	 }
      last15=j; 		   // save l value
      x15=x;			   // last x value
      }
   if ((j+49)==((j+49)/54)*54) {   // check if 54 divides l+49
      if ((x16>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last16-save16)/54; // difference in l values
	 save16=last16; 	   // save old l value
	 if (iter16!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x16)
	       max0[delta]=x16;
	    if (min0[delta]>x16)
	       min0[delta]=x16;
	    }
	 iter16=iter16+1;	   // increment period count
	 }
      last16=j; 		   // save l value
      x16=x;			   // last x value
      }
   if ((j+53)==((j+53)/54)*54) {   // check if 54 divides l+53
      if ((x17>0.0)&&(x<0.0)) {    // check for new period
	 delta=(last17-save17)/54; // difference in l values
	 save17=last17; 	   // save old l value
	 if (iter17!=0) {
	    if (delta<400)	   // histogram differences
	       hist0[delta]+=1;
	    else {
	       error[0]=3;
	       error[1]=delta;
	       goto zskip;
	       }
	    if (max0[delta]<x17)
	       max0[delta]=x17;
	    if (min0[delta]>x17)
	       min0[delta]=x17;
	    }
	 iter17=iter17+1;	   // increment period count
	 }
      last17=j; 		   // save l value
      x17=x;			   // last x value
      }
bskip:
   j=0;
   }
outx[0]=maxx;
outx[1]=minx;
outp[0]=maxp;
outp[1]=minp;
for (i=0; i<400; i++)
   maxi[i]=(unsigned int)(max0[i]/65536.0);
for (i=0; i<198140; i++)
   XG[i]=X[i];
zskip:
return(0);
}