~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

Linux Cross Reference
Tina6/tina-libs/tina/image/imgPrc_shade.c

Version: ~
Architecture: ~ [ i386 ] ~ [ alpha ] ~ [ m68k ] ~ [ mips ] ~ [ ppc ] ~ [ sparc ] ~ [ sparc64 ] ~

  1 /**********
  2  *
  3  * Copyright (c) 2003, Division of Imaging Science and Biomedical Engineering,
  4  * University of Manchester, UK.  All rights reserved.
  5  * 
  6  * Redistribution and use in source and binary forms, with or without modification, 
  7  * are permitted provided that the following conditions are met:
  8  * 
  9  *   . Redistributions of source code must retain the above copyright notice, 
 10  *     this list of conditions and the following disclaimer.
 11  *    
 12  *   . Redistributions in binary form must reproduce the above copyright notice,
 13  *     this list of conditions and the following disclaimer in the documentation 
 14  *     and/or other materials provided with the distribution.
 15  * 
 16  *   . Neither the name of the University of Manchester nor the names of its
 17  *     contributors may be used to endorse or promote products derived from this 
 18  *     software without specific prior written permission.
 19  * 
 20  * 
 21  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
 22  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 
 23  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 
 24  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 
 25  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
 26  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
 27  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 28  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 
 29  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 
 30  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 
 31  * POSSIBILITY OF SUCH DAMAGE.
 32  *
 33  **********
 34  *
 35  * Program :    TINA
 36  * File    :  $Source: /home/tina/cvs/tina-libs/tina/image/imgPrc_shade.c,v $
 37  * Date    :  $Date: 2009/03/19 18:36:26 $
 38  * Version :  $Revision: 1.7 $
 39  * CVS Id  :  $Id: imgPrc_shade.c,v 1.7 2009/03/19 18:36:26 paul Exp $
 40  *
 41  * Author  : Legacy TINA
 42  */
 43 
 44 /** 
 45  *  @file
 46  *  @brief Routines for applying a simple Lambertian illumination model to
 47  *  the generation of shaded images of a depth map. Includes routines which
 48  *  approximate the shading process as a convolution, with equivalent convolution kernel 
 49  *  generated in the Fourier domain. This leads on to attepts at shape from shading
 50  *  based upon deconvolution, which are largely unstable and inaccurate.
 51  *
 52  */
 53 
 54 #include "imgPrc_shade.h"
 55 
 56 #if HAVE_CONFIG_H
 57 #include <config.h>
 58 #endif
 59 
 60 #include <math.h>
 61 #include <tina/sys/sysDef.h>
 62 #include <tina/sys/sysPro.h>
 63 #include <tina/math/mathDef.h>
 64 #include <tina/math/mathPro.h>
 65 #include <tina/image/img_GenDef.h>
 66 #include <tina/image/img_GenPro.h>
 67 #include <tina/image/imgPrc_apply.h>
 68 #include <tina/image/imgPrc_grad.h>
 69 
 70 
 71 Imrect         *im_shading(Imrect * im, double slant, double tilt, double scale)
 72 {
 73         /* Shade an image */
 74 
 75         Imrect         *fim, *hgim, *vgim;
 76         Imregion       *roi;
 77         float          *row1, *row2, *row3;
 78         int             i, j;
 79         int             lx, ux, ly, uy;
 80         double          ctss, stss, cs;
 81 
 82         fim = im_re(im);
 83         if (fim == NULL)
 84                 return (NULL);
 85 
 86         roi = fim->region;
 87         if (roi == NULL)
 88                 return (NULL);
 89         lx = roi->lx;
 90         ux = roi->ux;
 91         ly = roi->ly;
 92         uy = roi->uy;
 93 
 94         hgim = imf_grad_h(fim);
 95         vgim = imf_grad_v(fim);
 96         ctss = cos(tilt) * sin(slant);
 97         stss = sin(tilt) * sin(slant);
 98         cs = cos(slant);
 99 
100         row1 = fvector_alloc(lx, ux);
101         row2 = fvector_alloc(lx, ux);
102         row3 = fvector_alloc(lx, ux);
103 
104         for (i = ly; i < uy; ++i)
105         {
106                 im_get_rowf(row1, hgim, i, lx, ux);
107                 im_get_rowf(row2, vgim, i, lx, ux);
108                 for (j = lx; j < ux; ++j)
109                 {
110                         row3[j] = (float) (scale * (row1[j] * ctss + row2[j] * stss + cs)
111                         / sqrt(row1[j] * row1[j] + row2[j] * row2[j] + 1.0));
112 
113                 }
114                 im_put_rowf(row3, fim, i, lx, ux);
115         }
116 
117         fvector_free(row1, lx);
118         fvector_free(row2, lx);
119         fvector_free(row3, lx);
120         im_free(hgim);
121         im_free(vgim);
122         return (fim);
123 }
124 
125 Imrect         *shade_conv(double slant, double tilt, Imregion * roi)
126 {
127         /* Compute effective convolution kernal in the fourier domain */
128 
129         Imrect         *im;
130         Complex        *row;
131         Complex         cgrad;
132         double          amplitude, phase;
133         int             lx, ux, ly, uy;
134         double          rangex, rangey;
135         int             i, j;
136 
137         if (roi == NULL)
138                 return (NULL);
139         lx = roi->lx;
140         ux = roi->ux;
141         ly = roi->ly;
142         uy = roi->uy;
143         rangex = roi->ux - roi->lx;
144         rangey = roi->uy - roi->ly;
145 
146         /* Added to prevent compiler warning PAB 04/12/2008 */
147         cgrad.ts_id=Complex_id;
148 
149         im = im_alloc(uy - ly, ux - lx, roi, complex_v);
150         row = zvector_alloc(lx, ux);
151 
152         for (i = ly; i < uy; ++i)
153         {
154                 for (j = lx; j < ux; ++j)
155                 {
156                         phase = -PIBY2;
157                         amplitude = sin(slant) * (sin(TWOPI * (j - lx) / rangex) * cos(tilt) +
158                                 sin(TWOPI * (i - ly) / rangey) * sin(tilt));
159                         /*
160                                     cgrad.y = amplitude * sin(phase);
161                                     cgrad.x = amplitude * cos(phase);
162                         */
163                         cgrad.y = -amplitude;
164                         cgrad.x = 0.0;
165                         row[j] = cgrad;
166                 }
167                 im_put_rowz(row, im, i, lx, ux);
168         }
169 
170         zvector_free(row, lx);
171         return (im);
172 }
173 
174 Imrect         *imz_fshade(Imrect * im1, double slant, double tilt)
175 /*
176  * compute the effective shading convolution for shallow surfaces in the
177  * fourier domain
178  */
179 {
180         Imrect         *im2;
181         Imregion       *roi;
182         Complex        *row1;
183         Complex         cgrad;
184         Complex        *row2;
185         double          amplitude, phase;
186         int             lx, ux, ly, uy;
187         double          rangex, rangey;
188         int             i, j;
189 
190         if (im1 == NULL || im1->vtype != complex_v)
191                 return (NULL);
192 
193         roi = im1->region;
194         if (roi == NULL)
195                 return (NULL);
196         lx = roi->lx;
197         ux = roi->ux;
198         ly = roi->ly;
199         uy = roi->uy;
200         rangex = roi->ux - roi->lx;
201         rangey = roi->uy - roi->ly;
202 
203         im2 = im_alloc(im1->height, im1->width, roi, complex_v);
204         row1 = zvector_alloc(lx, ux);
205         row2 = zvector_alloc(lx, ux);
206 
207         /* Added to prevent compiler warning PAB 04/12/2008 */
208         cgrad.ts_id=Complex_id;
209 
210         for (i = ly; i < uy; ++i)
211         {
212                 im_get_rowz(row1, im1, i, lx, ux);
213                 for (j = lx; j < ux; ++j)
214                 {
215                         cgrad.x = (row1[j]).x;
216                         cgrad.y = (row1[j]).y;
217                         if (cgrad.x * cgrad.x + cgrad.y * cgrad.y == 0.0)
218                                 amplitude = 0.0;
219                         else
220                                 amplitude = sqrt(cgrad.x * cgrad.x + cgrad.y * cgrad.y);
221 
222                         if (cgrad.x != 0.0)
223                                 phase = atan2(cgrad.y, cgrad.x);
224                         else if (cgrad.y > 0.0)
225                                 phase = PIBY2;
226                         else
227                                 phase = -PIBY2;
228                         phase -= PIBY2;
229                         amplitude *= sin(slant) * (sin(TWOPI * (j - lx) / rangex) * cos(tilt) +
230                                 sin(TWOPI * (i - ly) / rangey) * sin(tilt));
231                         cgrad.y = amplitude * sin(phase);
232                         cgrad.x = amplitude * cos(phase);
233                         row2[j] = cgrad;
234                 }
235                 im_put_rowz(row2, im2, i, lx, ux);
236         }
237 
238         zvector_free(row1, lx);
239         zvector_free(row2, lx);
240         return (im2);
241 }
242 
243 Imrect         *imz_fshape(Imrect * im1, double slant, double tilt, double limit)
244 {
245         Imrect         *im2;
246         Imregion       *roi;
247         Complex        *row1;
248         Complex         cgrad;
249         Complex        *row2;
250         double          amplitude, phase, fac;
251         int             lx, ux, ly, uy;
252         double          rangex, rangey;
253         double          ss, ct, st;
254         int             i, j;
255 
256         if (im1 == NULL || im1->vtype != complex_v)
257                 return (NULL);
258 
259         roi = im1->region;
260         if (roi == NULL)
261                 return (NULL);
262         lx = roi->lx;
263         ux = roi->ux;
264         ly = roi->ly;
265         uy = roi->uy;
266         rangex = roi->ux - roi->lx;
267         rangey = roi->uy - roi->ly;
268 
269         im2 = im_alloc(im1->height, im1->width, roi, complex_v);
270         row1 = zvector_alloc(lx, ux);
271         row2 = zvector_alloc(lx, ux);
272         ss = sin(slant);
273         ct = cos(tilt);
274         st = sin(tilt);
275 
276         /* Added to prevent compiler warning PAB 04/12/2008 */
277         cgrad.ts_id=Complex_id;
278 
279         for (i = ly; i < uy; ++i)
280         {
281                 im_get_rowz(row1, im1, i, lx, ux);
282                 for (j = lx; j < ux; ++j)
283                 {
284                         cgrad.x = (row1[j]).x;
285                         cgrad.y = (row1[j]).y;
286                         if (cgrad.x * cgrad.x + cgrad.y * cgrad.y == 0.0)
287                                 amplitude = 0.0;
288                         else
289                                 amplitude = sqrt(cgrad.x * cgrad.x + cgrad.y * cgrad.y);
290 
291                         if (cgrad.x != 0.0)
292                                 phase = atan2(cgrad.y, cgrad.x);
293                         else if (cgrad.y > 0.0)
294                                 phase = PIBY2;
295                         else
296                                 phase = -PIBY2;
297 
298                         phase += PIBY2;
299                         fac = ss * (sin(TWOPI * (j - lx) / rangex) * ct +
300                                     sin(TWOPI * (i - ly) / rangey) * st);
301                         if (fac > 0.0)
302                         {
303                                 if (fac < limit * ss)
304                                         fac = limit * ss;
305                         } else
306                         {
307                                 if (fac > -limit * ss)
308                                         fac = -limit * ss;
309                         }
310                         if (!(fabs(fac) > 0.00001))
311                                 format("bollocks\n");
312 
313                         amplitude /= fac;
314                         cgrad.y = amplitude * sin(phase);
315                         cgrad.x = amplitude * cos(phase);
316                         row2[j] = cgrad;
317                 }
318                 im_put_rowz(row2, im2, i, lx, ux);
319         }
320 
321         zvector_free(row1, lx);
322         zvector_free(row2, lx);
323         return (im2);
324 }
325 
326 Imrect         *imz_fxgrad(Imrect * im1)
327 /* compute x gradient of image in the fourier domain */
328 {
329         Imrect         *im2;
330         Imregion       *roi;
331         Complex        *row1;
332         Complex         cgrad;
333         Complex        *row2;
334         double          amplitude, phase;
335         int             lx, ux, ly, uy;
336         double          range;
337         int             i, j;
338 
339         if (im1 == NULL || im1->vtype != complex_v)
340                 return (NULL);
341 
342         roi = im1->region;
343         if (roi == NULL)
344                 return (NULL);
345         lx = roi->lx;
346         ux = roi->ux;
347         ly = roi->ly;
348         uy = roi->uy;
349         range = roi->ux - roi->lx;
350 
351         im2 = im_alloc(im1->height, im1->width, roi, complex_v);
352         row1 = zvector_alloc(lx, ux);
353         row2 = zvector_alloc(lx, ux);
354 
355         /* Added to prevent compiler warning PAB 04/12/2008 */
356         cgrad.ts_id=Complex_id;
357 
358         for (i = ly; i < uy; ++i)
359         {
360                 im_get_rowz(row1, im1, i, lx, ux);
361                 for (j = lx; j < ux; ++j)
362                 {
363                         cgrad.x = (row1[j]).x;
364                         cgrad.y = (row1[j]).y;
365                         if (cgrad.x * cgrad.x + cgrad.y * cgrad.y == 0.0)
366                                 amplitude = 0.0;
367                         else
368                                 amplitude = sqrt(cgrad.x * cgrad.x + cgrad.y * cgrad.y);
369 
370                         phase = atan2(cgrad.y, cgrad.x);
371                         phase -= PIBY2;
372                         amplitude *= sin(TWOPI * (j - lx) / range);
373                         cgrad.y = amplitude * sin(phase);
374                         cgrad.x = amplitude * cos(phase);
375                         row2[j] = cgrad;
376                 }
377                 im_put_rowz(row2, im2, i, lx, ux);
378         }
379 
380         zvector_free(row1, lx);
381         zvector_free(row2, lx);
382         return (im2);
383 }
384 
385 Imrect         *imz_fygrad(Imrect * im1)
386 /* compute y gradient of image in the fourier domain */
387 {
388         Imrect         *im2;
389         Imregion       *roi;
390         Complex        *row1;
391         Complex         cgrad;
392         Complex        *row2;
393         double          amplitude, phase;
394         int             lx, ux, ly, uy;
395         double          range;
396         int             i, j;
397 
398         if (im1 == NULL || im1->vtype != complex_v)
399                 return (NULL);
400 
401         roi = im1->region;
402         if (roi == NULL)
403                 return (NULL);
404         lx = roi->lx;
405         ux = roi->ux;
406         ly = roi->ly;
407         uy = roi->uy;
408         range = roi->uy - roi->ly;
409 
410         im2 = im_alloc(im1->height, im1->width, roi, complex_v);
411         row1 = zvector_alloc(lx, ux);
412         row2 = zvector_alloc(lx, ux);
413 
414         /* Added to prevent compiler warning PAB 04/12/2008 */
415         cgrad.ts_id=Complex_id;
416 
417         for (i = ly; i < uy; ++i)
418         {
419                 im_get_rowz(row1, im1, i, lx, ux);
420                 for (j = lx; j < ux; ++j)
421                 {
422                         cgrad.x = (row1[j]).x;
423                         cgrad.y = (row1[j]).y;
424                         if (cgrad.x * cgrad.x + cgrad.y * cgrad.y == 0.0)
425                                 amplitude = 0.0;
426                         else
427                                 amplitude = sqrt(cgrad.x * cgrad.x + cgrad.y * cgrad.y);
428 
429                         phase = atan2(cgrad.y, cgrad.x);
430                         phase -= PIBY2;
431                         amplitude *= sin(TWOPI * (i - ly) / range);
432                         cgrad.y = amplitude * sin(phase);
433                         cgrad.x = amplitude * cos(phase);
434                         row2[j] = cgrad;
435                 }
436                 im_put_rowz(row2, im2, i, lx, ux);
437         }
438 
439         zvector_free(row1, lx);
440         zvector_free(row2, lx);
441         return (im2);
442 }
443 

~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

This page was automatically generated by the LXR engine.
Visit the LXR main site for more information.