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_vcase1.cpp $
5 $LastChangedRevision: 85 $
6 $Date: 2014-07-02 16:07:08 -0400 (Wed, 02 Jul 2014) $
8 ***********************************************************
10 ------------------------------------------------------------------
12 VCASE1: VARIABLE SPACING CASE 1
13 PURPOSE: FOR VARIABLE SPACING, ICASE=1 (THE TWO LINES LIE BETWEEN
14 -90 AND -45 DEGREE LINES) FIND ALL POINTS THAT LIE
15 BETWEEN THE TWO LINES AND INTERPOLATE
17 ------------------------------------------------------------------
19 // bug79 - problem with fourier algorithm
29 #define cot(z) (cos(z)/sin(z));
31 void foru_class::vcase1(REAL ang1, REAL ang2, REAL* pj1, REAL* pj2)
68 uleft = (REAL) cot(ang1);
69 uright = (REAL) cot(ang2);
70 kh = Fourie.nsize2 / 2 - 1;
72 for (j = 1; j <= kh; j++)
74 // bug79 - swr 2/11/05
79 indl = (INTEGER) ceil(aleft - Fourie.epsln); // bug82 - changed to ceil - swr - 2/20/05
80 indh = (INTEGER) floor(aright + Fourie.epsln); // bug82 - changed to floor - swr - 2/20/05
83 a = j * Fourie.delp / Fourie.deld;
85 indmin = 2 * low; // bug79 - swr 2/11/05
87 distan = aright - aleft;
89 switch (Fourie.interp)
92 // INTERP = 1 NEAREST NEIGHBOR IN CARTESIAN SENSE
93 // INTERP = 4 F=(F1/D1+F2/D2+F3/D3+F4/D4)/(1/D1+1/D2+1/D3+1/D4)
99 yy[2] = (REAL) -low - 1;
100 yy[4] = (REAL) -low - 1;
101 xx[1] = yy[1] * uleft;
102 xx[2] = yy[2] * uleft;
103 xx[3] = yy[3] * uright;
104 xx[4] = yy[4] * uright;
106 for (ix = indl; ix <= indh; ix++)
107 { // bug79 - swr 2/11/05
108 xg = ix * Fourie.delp / Fourie.deld;
109 vint14(xx, yy, xg, yg, pj1, pj2, indmin, indmax, &gr, &gi); // changed to pass address of "gr" and "gi". Lajos, Feb 10, 2005
110 store(ix, iy, gr, gi);
113 break; // bug79 - swr 2/15/05
115 // INTERP = 2 NEAREST NEIGHBOR IN RADIAL SENSE
118 amid = (aleft + aright) / (REAL) 2.0;
119 ind = 2 * (lround(a)) + 1; // changed "(INTEGER)" to "lround" to match Fortran code. Lajos, Jan 20, 2005
120 for (ix = indl; ix <= indh; ix++)
121 { // bug79 - swr 2/11/05
124 { // bug79 - swr 2/15/05
133 store(ix, iy, gr, gi);
137 // INTERP = 3 LINEAR INTERPOLATION
140 if (t <= Consts.zero)
142 if ((1. - t) <= Consts.zero)
144 pr = pj1[indmin] * ((REAL) 1.0 - t) + pj1[indmax] * t;
145 px = pj1[indmin + 1] * ((REAL) 1.0 - t) + pj1[indmax + 1] * t;
146 qr = pj2[indmin] * ((REAL) 1.0 - t) + pj2[indmax] * t;
147 qi = pj2[indmin + 1] * ((REAL) 1.0 - t) + pj2[indmax + 1] * t;
149 for (ix = indl; ix <= indh; ix++)
150 { // bug79 - swr 2/11/05
151 s = (ix - aleft) / distan;
152 if (s <= Consts.zero)
154 if ((1.0 - s) <= Consts.zero)
156 gr = pr * ((REAL) 1.0 - s) + qr * s;
157 gi = px * ((REAL) 1.0 - s) + qi * s;
158 store(ix, iy, gr, gi);