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_rank.c,v $
37 * Date : $Date: 2008/12/04 15:28:02 $
38 * Version : $Revision: 1.7 $
39 * CVS Id : $Id: imgPrc_rank.c,v 1.7 2008/12/04 15:28:02 paul Exp $
40 *
41 * Author : Legacy TINA
42 */
43
44 /**
45 * @file
46 * @brief Rank filter in user defined square window.
47 *
48 *
49 * imi_rank: Generates an images where each pixel is the rank in a
50 * window of *size 2*size +1 with a real valued rank determined according to
51 * the provided noise estimate.
52 *
53 */
54
55 #include "imgPrc_rank.h"
56
57 #if HAVE_CONFIG_H
58 #include <config.h>
59 #endif
60
61 #include <math.h>
62 #include <tina/sys/sysDef.h>
63 #include <tina/sys/sysPro.h>
64 #include <tina/math/mathDef.h>
65 #include <tina/math/mathPro.h>
66 #include <tina/image/img_GenDef.h>
67 #include <tina/image/img_GenPro.h>
68
69
70 #define MAXRANGE 10 /* maximum of unrolled for loop NAT 3/9/2000 */
71 #define SOFT_RANK( x , y , z) (( (y)-(x) )< (z) ) ? (((y)-(x))/(z)) : 1.0
72
73
74 Imrect *imi_rank(Imrect * im, int size, double noise)
75 {
76 int lx, ux, ly, uy, i, j;
77 int m, icen, jcen;
78 int isize;
79 Imrect *fim;
80 int *rowpix, *newpix;
81 int *otherpix;
82 double cenpix;
83 double level;
84
85 if (im == NULL)
86 return (NULL);
87 lx = im->region->lx;
88 ux = im->region->ux;
89 ly = im->region->ly;
90 uy = im->region->uy;
91
92 fim = im_alloc(im->height, im->width, im->region, int_v);
93 for (i = ly; i < uy; i++)
94 {
95 newpix = IM_ROW(fim, i);
96 icen = i;
97 if (i < ly + size)
98 icen = ly + size;
99 if (i > uy - size - 1)
100 icen = uy - size - 1;
101 for (j = lx; j < ux; j++)
102 {
103 jcen = j;
104 if (j < lx + size)
105 jcen = lx + size;
106 if (j > ux - size - 1)
107 jcen = ux - size - 1;
108 cenpix = IM_INT(im, i, j) + 0.5*noise;
109 level = 0.0;
110 isize = icen + size + 1;
111 for (m = icen - size; m < isize; m++)
112 {
113 rowpix = IM_ROW(im, m);
114 otherpix = &rowpix[jcen - size - 1];
115 switch (size)
116 {
117 case 10:
118 if (*(++otherpix) < cenpix)
119 level += SOFT_RANK(*otherpix, cenpix, noise);
120 if (*(++otherpix) < cenpix)
121 level += SOFT_RANK(*otherpix, cenpix, noise);
122 case 9:
123 if (*(++otherpix) < cenpix)
124 level += SOFT_RANK(*otherpix, cenpix, noise);
125 if (*(++otherpix) < cenpix)
126 level += SOFT_RANK(*otherpix, cenpix, noise);
127 case 8:
128 if (*(++otherpix) < cenpix)
129 level += SOFT_RANK(*otherpix, cenpix, noise);
130 if (*(++otherpix) < cenpix)
131 level += SOFT_RANK(*otherpix, cenpix, noise);
132 case 7:
133 if (*(++otherpix) < cenpix)
134 level += SOFT_RANK(*otherpix, cenpix, noise);
135 if (*(++otherpix) < cenpix)
136 level += SOFT_RANK(*otherpix, cenpix, noise);
137 case 6:
138 if (*(++otherpix) < cenpix)
139 level += SOFT_RANK(*otherpix, cenpix, noise);
140 if (*(++otherpix) < cenpix)
141 level += SOFT_RANK(*otherpix, cenpix, noise);
142 case 5:
143 if (*(++otherpix) < cenpix)
144 level += SOFT_RANK(*otherpix, cenpix, noise);
145 if (*(++otherpix) < cenpix)
146 level += SOFT_RANK(*otherpix, cenpix, noise);
147 case 4:
148 if (*(++otherpix) < cenpix)
149 level += SOFT_RANK(*otherpix, cenpix, noise);
150 if (*(++otherpix) < cenpix)
151 level += SOFT_RANK(*otherpix, cenpix, noise);
152 case 3:
153 if (*(++otherpix) < cenpix)
154 level += SOFT_RANK(*otherpix, cenpix, noise);
155 if (*(++otherpix) < cenpix)
156 level += SOFT_RANK(*otherpix, cenpix, noise);
157 case 2:
158 if (*(++otherpix) < cenpix)
159 level += SOFT_RANK(*otherpix, cenpix, noise);
160 if (*(++otherpix) < cenpix)
161 level += SOFT_RANK(*otherpix, cenpix, noise);
162 case 1:
163 if (*(++otherpix) < cenpix)
164 level += SOFT_RANK(*otherpix, cenpix, noise);
165 if (*(++otherpix) < cenpix)
166 level += SOFT_RANK(*otherpix, cenpix, noise);
167 if (*(++otherpix) < cenpix)
168 level += SOFT_RANK(*otherpix, cenpix, noise);
169 }
170 }
171 newpix[j] = level;
172 }
173 }
174 return (fim);
175 }
176
177 /* For test purposes only NAT 3/9/2000 */
178 /* Commented out to prevent compiler warning PAB 04/12/2008 */
179 /*
180 static double soft_rank(float otherpix, double cenpix, double noise)
181 {
182 if (cenpix - otherpix < noise)
183 return ((cenpix - otherpix) / noise);
184 else
185 return (1.0);
186
187 }
188 */
189
190 Imrect *imf_rank(Imrect * im, int size, double noise)
191 {
192 int lx, ux, ly, uy, i, j;
193 int m, icen, jcen;
194 int isize;
195 Imrect *fim;
196 float *rowpix, *newpix;
197 float *otherpix;
198 double cenpix;
199 double level;
200
201 if (im == NULL)
202 return (NULL);
203 /*
204 * try to introduce a sensible radially symetric filter static float
205 * **weight=NULL; radial gaussian weighting with sigma = size/2.0
206 * significantly increases computation and required region size but
207 * does not improve final result. NAT 8/6/95 if (weight == NULL) {
208 * weight = farray_alloc(-size,-size,size+1,size+1); for
209 * (i=-size;i<size+1;i++) { for (j=-size;j<size+1;j++) { weight[i][j]
210 * = exp(-2.0*(i*i + j*j)/(double)(size*size)); } } }
211 */
212
213 lx = im->region->lx;
214 ux = im->region->ux;
215 ly = im->region->ly;
216 uy = im->region->uy;
217
218 fim = im_alloc(im->height, im->width, im->region, float_v);
219 for (i = ly; i < uy; i++)
220 {
221 newpix = IM_ROW(fim, i);
222 icen = i;
223 if (i < ly + size)
224 icen = ly + size;
225 if (i > uy - size - 1)
226 icen = uy - size - 1;
227 for (j = lx; j < ux; j++)
228 {
229 jcen = j;
230 if (j < lx + size)
231 jcen = lx + size;
232 if (j > ux - size - 1)
233 jcen = ux - size - 1;
234 cenpix = IM_FLOAT(im, i, j) + 0.5*noise;
235 level = 0.0;
236 isize = icen + size + 1;
237
238 for (m = icen - size; m < isize; m++)
239 {
240 rowpix = IM_ROW(im, m);
241 otherpix = &rowpix[jcen - size - 1];
242 switch (size) /* unrolled loop for speed */
243 {
244 case 10:
245 if (*(++otherpix) < cenpix )
246 level += SOFT_RANK(*otherpix, cenpix, noise);
247 if (*(++otherpix) < cenpix )
248 level += SOFT_RANK(*otherpix, cenpix, noise);
249 case 9:
250 if (*(++otherpix) < cenpix )
251 level += SOFT_RANK(*otherpix, cenpix, noise);
252 if (*(++otherpix) < cenpix )
253 level += SOFT_RANK(*otherpix, cenpix, noise);
254 case 8:
255 if (*(++otherpix) < cenpix )
256 level += SOFT_RANK(*otherpix, cenpix, noise);
257 if (*(++otherpix) < cenpix )
258 level += SOFT_RANK(*otherpix, cenpix, noise);
259 case 7:
260 if (*(++otherpix) < cenpix )
261 level += SOFT_RANK(*otherpix, cenpix, noise);
262 if (*(++otherpix) < cenpix )
263 level += SOFT_RANK(*otherpix, cenpix, noise);
264 case 5:
265 if (*(++otherpix) < cenpix )
266 level += SOFT_RANK(*otherpix, cenpix, noise);
267 if (*(++otherpix) < cenpix )
268 level += SOFT_RANK(*otherpix, cenpix, noise);
269 case 4:
270 if (*(++otherpix) < cenpix )
271 level += SOFT_RANK(*otherpix, cenpix, noise);
272 if (*(++otherpix) < cenpix )
273 level += SOFT_RANK(*otherpix, cenpix, noise);
274 case 3:
275 if (*(++otherpix) < cenpix )
276 level += SOFT_RANK(*otherpix, cenpix, noise);
277 if (*(++otherpix) < cenpix )
278 level += SOFT_RANK(*otherpix, cenpix, noise);
279 case 2:
280 if (*(++otherpix) < cenpix )
281 level += SOFT_RANK(*otherpix, cenpix, noise);
282 if (*(++otherpix) < cenpix )
283 level += SOFT_RANK(*otherpix, cenpix, noise);
284 case 1:
285 if (*(++otherpix) < cenpix )
286 level += SOFT_RANK(*otherpix, cenpix, noise);
287 if (*(++otherpix) < cenpix )
288 level += SOFT_RANK(*otherpix, cenpix, noise);
289 if (*(++otherpix) < cenpix )
290 level += SOFT_RANK(*otherpix, cenpix, noise);
291 }
292 }
293 newpix[j] = level;
294 }
295 }
296 return (fim);
297 }
298
299 Imrect *im_rank(Imrect * im, int range, double noise)
300 {
301 Imrect *rim;
302 if (im == NULL || range > MAXRANGE)
303 return (NULL);
304 switch (im->vtype)
305 {
306 case uchar_v:
307 case char_v:
308 case short_v:
309 case ushort_v:
310 im = im_cast(im, int_v);
311 rim = imi_rank(im, range, noise);
312 im_free(im);
313 return (rim);
314 case int_v:
315 return (imi_rank(im, range, noise));
316 case double_v:
317 im = im_cast(im, float_v);
318 rim = imf_rank(im, range, noise);
319 im_free(im);
320 return (rim);
321 case float_v:
322 return (imf_rank(im, range, noise));
323 case complex_v:
324 return (NULL);
325 default:
326 return (NULL);
327 }
328 }
329
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.