1 #include <math.h>
2 #include <tina/sys.h>
3 #include <tina/sysfuncs.h>
4 #include <tina/math.h>
5 #include <tina/mathfuncs.h>
6
7 Imrect *im_shading(Imrect *im,double slant,double tilt,double scale)
8 /* shade an image */
9 {
10 Imrect *fim,*hgim,*vgim,*im_re();
11 Imrect *imf_grad_h(),*imf_grad_v();
12 Imregion *roi;
13 float *row1,*row2,*row3;
14 int i,j;
15 int lx,ux,ly,uy;
16 double ctss,stss,cs;
17
18 fim = im_re(im);
19 if (fim == NULL)
20 return (NULL);
21
22 roi = fim->region;
23 if (roi == NULL)
24 return (NULL);
25 lx = roi->lx;
26 ux = roi->ux;
27 ly = roi->ly;
28 uy = roi->uy;
29
30 hgim = imf_grad_h(fim);
31 vgim = imf_grad_v(fim);
32 ctss = cos(tilt)*sin(slant);
33 stss = sin(tilt)*sin(slant);
34 cs = cos(slant);
35
36 row1 = fvector_alloc(lx, ux);
37 row2 = fvector_alloc(lx, ux);
38 row3 = fvector_alloc(lx, ux);
39
40 for (i = ly; i < uy; ++i)
41 {
42 im_get_rowf(row1, hgim, i, lx, ux);
43 im_get_rowf(row2, vgim, i, lx, ux);
44 for (j = lx; j < ux; ++j)
45 {
46 row3[j] = (float) (scale*(row1[j]*ctss + row2[j]*stss + cs)
47 / sqrt(row1[j]*row1[j] + row2[j]*row2[j] + 1.0));
48
49 }
50 im_put_rowf(row3, fim, i, lx, ux);
51 }
52
53 fvector_free(row1, lx);
54 fvector_free(row2, lx);
55 fvector_free(row3, lx);
56 im_free(hgim);
57 im_free(vgim);
58 return (fim);
59 }
60
61 Imrect *shade_conv(double slant,double tilt,Imregion *roi)
62 /* compute effective convolution kernal in the fourier domain */
63 {
64 Imrect *im;
65 Complex *row;
66 Complex cgrad;
67 double amplitude,phase;
68 int lx, ux, ly, uy;
69 double rangex,rangey;
70 int i, j;
71
72 if (roi == NULL)
73 return (NULL);
74 lx = roi->lx;
75 ux = roi->ux;
76 ly = roi->ly;
77 uy = roi->uy;
78 rangex = roi->ux-roi->lx;
79 rangey = roi->uy-roi->ly;
80
81 im = im_alloc(uy-ly, ux-lx, roi, complex_v);
82 row = zvector_alloc(lx, ux);
83
84 for (i = ly; i < uy; ++i)
85 {
86 for (j = lx; j < ux; ++j)
87 {
88 phase = -PIBY2;
89 amplitude = sin(slant)*(sin(TWOPI*(j-lx)/rangex)*cos(tilt) +
90 sin(TWOPI*(i-ly)/rangey)*sin(tilt));
91 /*
92 cgrad.y = amplitude * sin(phase);
93 cgrad.x = amplitude * cos(phase);
94 */
95 cgrad.y = -amplitude;
96 cgrad.x = 0.0;
97 row[j] = cgrad;
98 }
99 im_put_rowz(row, im, i, lx, ux);
100 }
101
102 zvector_free(row, lx);
103 return (im);
104 }
105
106 Imrect *imz_fshade(Imrect *im1,double slant,double tilt)
107 /* compute the effective shading convolution for shallow
108 surfaces in the fourier domain */
109 {
110 Imrect *im2;
111 Imregion *roi;
112 Complex *row1;
113 Complex cgrad;
114 Complex *row2;
115 double amplitude,phase;
116 int lx, ux, ly, uy;
117 double rangex,rangey;
118 int i, j;
119
120 if (im1 == NULL||im1->vtype!=complex_v)
121 return (NULL);
122
123 roi = im1->region;
124 if (roi == NULL)
125 return (NULL);
126 lx = roi->lx;
127 ux = roi->ux;
128 ly = roi->ly;
129 uy = roi->uy;
130 rangex = roi->ux-roi->lx;
131 rangey = roi->uy-roi->ly;
132
133 im2 = im_alloc(im1->height, im1->width, roi, complex_v);
134 row1 = zvector_alloc(lx, ux);
135 row2 = zvector_alloc(lx, ux);
136
137 for (i = ly; i < uy; ++i)
138 {
139 im_get_rowz(row1, im1, i, lx, ux);
140 for (j = lx; j < ux; ++j)
141 {
142 cgrad.x = (row1[j]).x;
143 cgrad.y = (row1[j]).y;
144 if (cgrad.x*cgrad.x + cgrad.y*cgrad.y==0.0)
145 amplitude = 0.0;
146 else
147 amplitude = sqrt(cgrad.x*cgrad.x + cgrad.y*cgrad.y);
148
149 if (cgrad.x!=0.0)
150 phase = atan2(cgrad.y,cgrad.x);
151 else if (cgrad.y > 0.0 ) phase = PIBY2;
152 else phase = -PIBY2;
153 phase -= PIBY2;
154 amplitude *= sin(slant)*(sin(TWOPI*(j-lx)/rangex)*cos(tilt) +
155 sin(TWOPI*(i-ly)/rangey)*sin(tilt));
156 cgrad.y = amplitude * sin(phase);
157 cgrad.x = amplitude * cos(phase);
158 row2[j] = cgrad;
159 }
160 im_put_rowz(row2, im2, i, lx, ux);
161 }
162
163 zvector_free(row1, lx);
164 zvector_free(row2, lx);
165 return (im2);
166 }
167
168 Imrect *imz_fshape(Imrect *im1,double slant,double tilt,double limit)
169 {
170 Imrect *im2;
171 Imregion *roi;
172 Complex *row1;
173 Complex cgrad;
174 Complex *row2;
175 double amplitude,phase,fac;
176 int lx, ux, ly, uy;
177 double rangex,rangey;
178 double ss,ct,st;
179 int i, j;
180
181 if (im1 == NULL||im1->vtype!=complex_v)
182 return (NULL);
183
184 roi = im1->region;
185 if (roi == NULL)
186 return (NULL);
187 lx = roi->lx;
188 ux = roi->ux;
189 ly = roi->ly;
190 uy = roi->uy;
191 rangex = roi->ux-roi->lx;
192 rangey = roi->uy-roi->ly;
193
194 im2 = im_alloc(im1->height, im1->width, roi, complex_v);
195 row1 = zvector_alloc(lx, ux);
196 row2 = zvector_alloc(lx, ux);
197 ss = sin(slant);
198 ct = cos(tilt);
199 st = sin(tilt);
200
201 for (i = ly; i < uy; ++i)
202 {
203 im_get_rowz(row1, im1, i, lx, ux);
204 for (j = lx; j < ux; ++j)
205 {
206 cgrad.x = (row1[j]).x;
207 cgrad.y = (row1[j]).y;
208 if (cgrad.x*cgrad.x + cgrad.y*cgrad.y==0.0)
209 amplitude = 0.0;
210 else
211 amplitude = sqrt(cgrad.x*cgrad.x + cgrad.y*cgrad.y);
212
213 if (cgrad.x!=0.0)
214 phase = atan2(cgrad.y,cgrad.x);
215 else if (cgrad.y > 0.0 ) phase = PIBY2;
216 else phase = -PIBY2;
217
218 phase += PIBY2;
219 fac = ss*(sin(TWOPI*(j-lx)/rangex)*ct +
220 sin(TWOPI*(i-ly)/rangey)*st);
221 if (fac>0.0)
222 {
223 if (fac<limit*ss) fac = limit*ss;
224 }
225 else
226 {
227 if (fac>-limit*ss) fac = -limit*ss;
228 }
229 if (!(fabs(fac)>0.00001))
230 format("bollocks\n");
231
232 amplitude /= fac;
233 cgrad.y = amplitude * sin(phase);
234 cgrad.x = amplitude * cos(phase);
235 row2[j] = cgrad;
236 }
237 im_put_rowz(row2, im2, i, lx, ux);
238 }
239
240 zvector_free(row1, lx);
241 zvector_free(row2, lx);
242 return (im2);
243 }
244
245 Imrect *imz_fxgrad(Imrect *im1)
246 /* compute x gradient of image in the fourier domain */
247 {
248 Imrect *im2;
249 Imregion *roi;
250 Complex *row1;
251 Complex cgrad;
252 Complex *row2;
253 double amplitude,phase;
254 int lx, ux, ly, uy;
255 double range;
256 int i, j;
257
258 if (im1 == NULL||im1->vtype!=complex_v)
259 return (NULL);
260
261 roi = im1->region;
262 if (roi == NULL)
263 return (NULL);
264 lx = roi->lx;
265 ux = roi->ux;
266 ly = roi->ly;
267 uy = roi->uy;
268 range = roi->ux-roi->lx;
269
270 im2 = im_alloc(im1->height, im1->width, roi, complex_v);
271 row1 = zvector_alloc(lx, ux);
272 row2 = zvector_alloc(lx, ux);
273
274 for (i = ly; i < uy; ++i)
275 {
276 im_get_rowz(row1, im1, i, lx, ux);
277 for (j = lx; j < ux; ++j)
278 {
279 cgrad.x = (row1[j]).x;
280 cgrad.y = (row1[j]).y;
281 if (cgrad.x*cgrad.x + cgrad.y*cgrad.y==0.0)
282 amplitude = 0.0;
283 else
284 amplitude = sqrt(cgrad.x*cgrad.x + cgrad.y*cgrad.y);
285
286 phase = atan2(cgrad.y,cgrad.x);
287 phase -= PIBY2;
288 amplitude *= sin(TWOPI*(j-lx)/range);
289 cgrad.y = amplitude * sin(phase);
290 cgrad.x = amplitude * cos(phase);
291 row2[j] = cgrad;
292 }
293 im_put_rowz(row2, im2, i, lx, ux);
294 }
295
296 zvector_free(row1, lx);
297 zvector_free(row2, lx);
298 return (im2);
299 }
300
301 Imrect *imz_fygrad(Imrect *im1)
302 /* compute y gradient of image in the fourier domain */
303 {
304 Imrect *im2;
305 Imregion *roi;
306 Complex *row1;
307 Complex cgrad;
308 Complex *row2;
309 double amplitude,phase;
310 int lx, ux, ly, uy;
311 double range;
312 int i, j;
313
314 if (im1 == NULL||im1->vtype!=complex_v)
315 return (NULL);
316
317 roi = im1->region;
318 if (roi == NULL)
319 return (NULL);
320 lx = roi->lx;
321 ux = roi->ux;
322 ly = roi->ly;
323 uy = roi->uy;
324 range = roi->uy-roi->ly;
325
326 im2 = im_alloc(im1->height, im1->width, roi, complex_v);
327 row1 = zvector_alloc(lx, ux);
328 row2 = zvector_alloc(lx, ux);
329
330 for (i = ly; i < uy; ++i)
331 {
332 im_get_rowz(row1, im1, i, lx, ux);
333 for (j = lx; j < ux; ++j)
334 {
335 cgrad.x = (row1[j]).x;
336 cgrad.y = (row1[j]).y;
337 if (cgrad.x*cgrad.x + cgrad.y*cgrad.y==0.0)
338 amplitude = 0.0;
339 else
340 amplitude = sqrt(cgrad.x*cgrad.x + cgrad.y*cgrad.y);
341
342 phase = atan2(cgrad.y,cgrad.x);
343 phase -= PIBY2;
344 amplitude *= sin(TWOPI*(i-ly)/range);
345 cgrad.y = amplitude * sin(phase);
346 cgrad.x = amplitude * cos(phase);
347 row2[j] = cgrad;
348 }
349 im_put_rowz(row2, im2, i, lx, ux);
350 }
351
352 zvector_free(row1, lx);
353 zvector_free(row2, lx);
354 return (im2);
355 }
356
357
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.