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

Linux Cross Reference
Tina6/tina-libs/tina/image/imgPrc_lsf.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_lsf.c,v $
 37  * Date    :  $Date: 2003/09/22 16:09:02 $
 38  * Version :  $Revision: 1.5 $
 39  * CVS Id  :  $Id: imgPrc_lsf.c,v 1.5 2003/09/22 16:09:02 tony Exp $
 40  *
 41  * Author  : Legacy TINA
 42  */
 43 
 44 /** 
 45  *  @file
 46  *  @brief Linear sequential filter functions for smoothing of images with large exponential
 47  *  kenels. Otherwise known as IIR and FIR filters.
 48  *
 49  */
 50 
 51 #include "imgPrc_lsf.h"
 52 
 53 #if HAVE_CONFIG_H
 54 #include <config.h>
 55 #endif
 56 
 57 #include <math.h>
 58 #include <tina/sys/sysDef.h>
 59 #include <tina/sys/sysPro.h>
 60 #include <tina/math/mathDef.h>
 61 #include <tina/math/mathPro.h>
 62 #include <tina/image/img_GenDef.h>
 63 #include <tina/image/img_GenPro.h>
 64 
 65 static void     lsf_smooth(int n1, int n2, float *u1, float *u2, double sigma, float *weight)
 66 {
 67         int             i;
 68         float          *v1 = fvector_alloc(n1, n2);
 69         float          *v2 = fvector_alloc(n1, n2);
 70         double          a, b = exp(-1.0 / sigma), ab;
 71 
 72         a = (1 - b) / (1 + b);  /** makes filter sum to 1 **/
 73         ab = a * b;
 74 
 75 #pragma omp parallel sections num_threads(2)
 76         {
 77             #pragma omp section
 78             {
 79                 v1[n1] = (float) (a * u1[n1]);
 80                 for (i = n1 + 1; i < n2; i++)
 81                     v1[i] = (float) (a * u1[i] + b * v1[i - 1]);
 82             }
 83 
 84             #pragma omp section
 85             {
 86                 v2[n2 - 1] = (float) 0.0;
 87                 for (i = n2 - 2; i >= n1; i--)
 88                     v2[i] = (float) (ab * u1[i + 1] + b * v2[i + 1]);
 89             }
 90         }
 91 
 92         for (i = n1; i < n2; i++)
 93                 u2[i] = v1[i] + v2[i];
 94 
 95         if (weight != NULL)
 96                 for (i = n1; i < n2; i++)
 97                         u2[i] /= weight[i];
 98 
 99         fvector_free(v1, n1);
100         fvector_free(v2, n1);
101 }
102 
103 static void     lsf_smooth_side(int n1, int n2, float *u1, float *u2, double sigma, float *weight, int side)
104 {
105         int             i;
106         double          a, b = exp(-1.0 / sigma);
107 
108         a = 1 - b;              /** makes filter sum to 1 **/
109 
110         if (side == 1)
111         {
112                 u2[n1] = (float) (a * u1[n1]);
113                 for (i = n1 + 1; i < n2; i++)
114                         u2[i] = (float) (a * u1[i] + b * u2[i - 1]);
115         } else if (side == -1)
116         {
117                 u2[n2 - 1] = (float) (a * u1[n2 - 1]);
118                 for (i = n2 - 2; i >= n1; i--)
119                         u2[i] = (float) (a * u1[i] + b * u2[i + 1]);
120         } else
121         {
122                 for (i = n1; i < n2; i++)
123                         u2[i] = (float) u1[i];
124         }
125 
126         if (weight != NULL)
127                 for (i = n1; i < n2; i++)
128                         u2[i] /= weight[i];
129 }
130 
131 Imrect         *imf_lsf_smooth_serial(Imrect * im1, double sigma)
132 {
133         Imrect         *im2;
134         Imrect         *im3;
135         Imregion       *roi;
136         float          *row1, *row2, *col1, *col2, *ones, *weight;
137         int             lx, ux, ly, uy;
138         int             i;
139 
140         if (im1 == NULL)
141                 return (NULL);
142 
143         roi = im1->region;
144         if (roi == NULL)
145                 return (NULL);
146         lx = roi->lx;
147         ux = roi->ux;
148         ly = roi->ly;
149         uy = roi->uy;
150 
151 
152         ones = fvector_alloc(lx, ux);
153         weight = fvector_alloc(lx, ux);
154         for (i = lx; i < ux; ++i)
155                 ones[i] = (float) 1.0;
156         lsf_smooth(lx, ux, ones, weight, sigma, (float *) NULL);
157         fvector_free(ones, lx);
158 
159         im2 = im_alloc(im1->height, im1->width, roi, float_v);
160         row1 = fvector_alloc(lx, ux);
161         row2 = fvector_alloc(lx, ux);
162 
163         for (i = ly; i < uy; ++i)
164         {
165                 im_get_rowf(row1, im1, i, lx, ux);
166                 lsf_smooth(lx, ux, row1, row2, sigma, weight);
167                 im_put_rowf(row2, im2, i, lx, ux);
168         }
169         fvector_free(row1, lx);
170         fvector_free(row2, lx);
171         fvector_free(weight, lx);
172 
173         ones = fvector_alloc(ly, uy);
174         weight = fvector_alloc(ly, uy);
175         for (i = ly; i < uy; ++i)
176                 ones[i] = (float) 1.0;
177         lsf_smooth(ly, uy, ones, weight, sigma, (float *) NULL);
178         fvector_free(ones, ly);
179 
180         im3 = im_alloc(im1->height, im1->width, roi, float_v);
181         col1 = fvector_alloc(ly, uy);
182         col2 = fvector_alloc(ly, uy);
183         for (i = lx; i < ux; ++i)
184         {
185                 im_get_colf(col1, im2, i, ly, uy);
186                 lsf_smooth(ly, uy, col1, col2, sigma, weight);
187                 im_put_colf(col2, im3, i, ly, uy);
188         }
189         fvector_free(col1, ly);
190         fvector_free(col2, ly);
191         fvector_free(weight, ly);
192 
193         im_free(im2);
194 
195         return (im3);
196 }
197 
198 
199 Imrect         *imf_lsf_smooth_parallel(Imrect * im1, double sigma)
200 {
201         Imrect         *im2;
202         Imrect         *im3;
203         Imregion       *roi;
204         float          *row1, *row2, *col1, *col2, *ones, *weight;
205         int             lx, ux, ly, uy;
206         int             i;
207 
208         if (im1 == NULL)
209                 return (NULL);
210 
211         roi = im1->region;
212         if (roi == NULL)
213                 return (NULL);
214         lx = roi->lx;
215         ux = roi->ux;
216         ly = roi->ly;
217         uy = roi->uy;
218 
219 
220         ones = fvector_alloc(lx, ux);
221         weight = fvector_alloc(lx, ux);
222         for (i = lx; i < ux; ++i)
223                 ones[i] = (float) 1.0;
224         lsf_smooth(lx, ux, ones, weight, sigma, (float *) NULL);
225         fvector_free(ones, lx);
226 
227         im2 = im_alloc(im1->height, im1->width, roi, float_v);
228 
229 #pragma omp parallel for private(row1, row2)
230         for (i = ly; i < uy; ++i)
231         {
232                 row1 = fvector_alloc(lx, ux);
233                 row2 = fvector_alloc(lx, ux);
234 
235                 im_get_rowf(row1, im1, i, lx, ux);
236                 lsf_smooth(lx, ux, row1, row2, sigma, weight);
237                 im_put_rowf(row2, im2, i, lx, ux);
238 
239                 fvector_free(row1, lx);
240                 fvector_free(row2, lx);
241         }
242 
243         fvector_free(weight, lx);
244 
245         ones = fvector_alloc(ly, uy);
246         weight = fvector_alloc(ly, uy);
247         for (i = ly; i < uy; ++i)
248                 ones[i] = (float) 1.0;
249         lsf_smooth(ly, uy, ones, weight, sigma, (float *) NULL);
250         fvector_free(ones, ly);
251 
252         im3 = im_alloc(im1->height, im1->width, roi, float_v);
253 
254 #pragma omp parallel for private(col1, col2)
255         for (i = lx; i < ux; ++i)
256         {
257                 col1 = fvector_alloc(ly, uy);
258                 col2 = fvector_alloc(ly, uy);
259  
260                 im_get_colf(col1, im2, i, ly, uy);
261                 lsf_smooth(ly, uy, col1, col2, sigma, weight);
262                 im_put_colf(col2, im3, i, ly, uy);
263 
264                 fvector_free(col1, ly);
265                 fvector_free(col2, ly);
266         }
267 
268         fvector_free(weight, ly);
269 
270         im_free(im2);
271 
272         return (im3);
273 }
274 
275 
276 Imrect         *imf_lsf_smooth(Imrect * im1, double sigma)
277 {
278     return imf_lsf_smooth_parallel(im1, sigma);
279 }
280 
281 
282 Imrect         *imf_lsf_smooth_quad(Imrect * im1, double sigma, int sidex, int sidey)
283 {
284         Imrect         *im2;
285         Imrect         *im3;
286         Imregion       *roi;
287         float          *row1, *row2, *col1, *col2, *ones, *weight;
288         int             lx, ux, ly, uy;
289         int             i;
290 
291         if (im1 == NULL)
292                 return (NULL);
293 
294         roi = im1->region;
295         if (roi == NULL)
296                 return (NULL);
297         lx = roi->lx;
298         ux = roi->ux;
299         ly = roi->ly;
300         uy = roi->uy;
301 
302 
303         ones = fvector_alloc(lx, ux);
304         weight = fvector_alloc(lx, ux);
305         for (i = lx; i < ux; ++i)
306                 ones[i] = (float) 1.0;
307         lsf_smooth_side(lx, ux, ones, weight, sigma, (float *) NULL, sidex);
308         fvector_free(ones, lx);
309 
310         im2 = im_alloc(im1->height, im1->width, roi, float_v);
311         row1 = fvector_alloc(lx, ux);
312         row2 = fvector_alloc(lx, ux);
313         for (i = ly; i < uy; ++i)
314         {
315                 im_get_rowf(row1, im1, i, lx, ux);
316                 lsf_smooth_side(lx, ux, row1, row2, sigma, weight, sidex);
317                 im_put_rowf(row2, im2, i, lx, ux);
318         }
319         fvector_free(row1, lx);
320         fvector_free(row2, lx);
321         fvector_free(weight, lx);
322 
323         ones = fvector_alloc(ly, uy);
324         weight = fvector_alloc(ly, uy);
325         for (i = ly; i < uy; ++i)
326                 ones[i] = (float) 1.0;
327         lsf_smooth_side(ly, uy, ones, weight, sigma, (float *) NULL, sidey);
328         fvector_free(ones, ly);
329 
330         im3 = im_alloc(im1->height, im1->width, roi, float_v);
331         col1 = fvector_alloc(ly, uy);
332         col2 = fvector_alloc(ly, uy);
333         for (i = lx; i < ux; ++i)
334         {
335                 im_get_colf(col1, im2, i, ly, uy);
336                 lsf_smooth_side(ly, uy, col1, col2, sigma, weight, sidey);
337                 im_put_colf(col2, im3, i, ly, uy);
338         }
339         fvector_free(col1, ly);
340         fvector_free(col2, ly);
341         fvector_free(weight, ly);
342 
343         im_free(im2);
344 
345         return (im3);
346 }
347 

~ [ 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.