2 ***********************************************************
3 $SNARK_Header: S N A R K 1 4 - A PICTURE RECONSTRUCTION PROGRAM $
4 $HeadURL: svn://dig.cs.gc.cuny.edu/snark/trunk/src/snark/foru_unif.cpp $
5 $LastChangedRevision: 85 $
6 $Date: 2014-07-02 16:07:08 -0400 (Wed, 02 Jul 2014) $
8 ***********************************************************
10 PURPOSE: FOR EACH CASE (ICASE=1,...,6) FIND ALL POINTS LIE BETWEEN
11 THE TWO RAYS AND CALL UINT TO CALAULATE INTERPOLATED
24 #define cot(z) (cos(z)/sin(z));
26 void foru_class::unif(INTEGER icase, REAL ang1, REAL ang2, REAL* nold, REAL* new1)
28 INTEGER kh, ix, iy, j, khd, kk;
29 REAL uleft, uright, uhigh, ulow, xhigh, xlow, yhigh, ylow;
31 kh = Fourie.nsize2 / 2 - 1;
38 uleft = (REAL) cot(ang1)
40 uright = (REAL) cot(ang2)
42 for (j = 1; j <= kh; j++)
45 xlow = uleft * iy - Fourie.epsln;
46 xhigh = uright * iy + Fourie.epsln;
47 ix = lround(xlow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
52 uint1(icase, ix, iy, ang1, ang2, nold, new1);
61 uleft = (REAL) cot(ang1)
63 for (j = 1; j <= kh; j++)
66 xlow = uleft * iy - Fourie.epsln;
67 ix = lround(xlow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
72 uint1(icase, ix, iy, ang1, ang2, nold, new1);
78 uhigh = (REAL) tan(ang2);
79 for (ix = 1; ix <= kh; ix++)
81 yhigh = ix * uhigh + Fourie.epsln;
85 uint1(icase, ix, iy, ang1, ang2, nold, new1);
93 ulow = (REAL) tan(ang1);
94 uhigh = (REAL) tan(ang2);
95 for (ix = 1; ix <= kh; ix++)
97 ylow = ulow * ix - Fourie.epsln;
98 yhigh = uhigh * ix + Fourie.epsln;
99 iy = lround(ylow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
104 uint1(icase, ix, iy, ang1, ang2, nold, new1);
114 ulow = (REAL) tan(ang1);
115 for (ix = 1; ix <= kh; ix++)
116 { // !!! ix = 0 does not work mk
117 ylow = ulow * ix - Fourie.epsln;
118 iy = lround(ylow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
120 { // changed ">=" to "<=" to match Fortran code. Lajos, Jan 13, 2005
123 uint1(icase, ix, iy, ang1, ang2, nold, new1);
129 uleft = (REAL) cot(ang2)
132 for (iy = 1; iy <= kh; iy++)
134 xlow = uleft * iy - Fourie.epsln;
138 uint1(icase, ix, iy, ang1, ang2, nold, new1);
148 uleft = (REAL) cot(ang2)
150 uright = (REAL) cot(ang1)
152 for (iy = 1; iy <= kh; iy++)
154 xlow = uleft * iy - Fourie.epsln;
155 xhigh = uright * iy + Fourie.epsln;
156 ix = lround(xlow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
161 uint1(icase, ix, iy, ang1, ang2, nold, new1);
170 // CONJUGATE THE SECOND PROJECTION
172 khd = Fourie.nsize1 + 2;
173 for (kk = 1; kk < khd; kk += 2)
175 new1[kk] = -new1[kk];
178 uleft = (REAL) cot(ang2)
180 uright = (REAL) cot(ang1)
182 for (iy = 1; iy <= kh; iy++)
184 xlow = uleft * iy - Fourie.epsln;
185 xhigh = uright * iy + Fourie.epsln;
186 ix = lround(xlow); // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
191 uint1(icase, ix, iy, ang1, ang2 + Consts.pi, nold, new1);
199 khd = Fourie.nsize1 + 2;
200 for (kk = 1; kk < khd; kk += 2)
202 new1[kk] = -new1[kk];