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

Linux Cross Reference
Tina6/tina-libs/tina/image/imgPrc_combine.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_combine.c,v $
 37  * Date    :  $Date: 2004/08/05 14:32:44 $
 38  * Version :  $Revision: 1.8 $
 39  * CVS Id  :  $Id: imgPrc_combine.c,v 1.8 2004/08/05 14:32:44 neil Exp $
 40  *
 41  * Author  : Legacy TINA
 42  */
 43 /** 
 44  *  @file
 45  *  @brief Simple arithmetic manipulation of two images and via a user specified per-pixel function.
 46  *
 47  * im_combine: combining two images, general case
 48  * all others: combining two floating point images, in place version, and special cases
 49  *
 50  */
 51 
 52 #include "imgPrc_combine.h"
 53 
 54 #if HAVE_CONFIG_H
 55 #include <config.h>
 56 #endif
 57 
 58 #include <math.h>
 59 #include <tina/sys/sysDef.h>
 60 #include <tina/sys/sysPro.h>
 61 #include <tina/math/mathDef.h>
 62 #include <tina/math/mathPro.h>
 63 #include <tina/image/img_GenDef.h>
 64 #include <tina/image/img_GenPro.h>
 65 #include <tina/image/imgPrc_suptype.h>
 66 
 67 
 68 /**
 69  * @brief Apply a user-specified per-pixel function to a pair of images to produce an integer image.
 70  * @param im1 Image to be combined.
 71  * @param im2 Image to be combined.
 72  * @param func Per-pixel combination function to use.
 73  * @param data Pointer to any additional data required by the combination function
 74  * @return im3 Integer result of the image combination
 75  *
 76  * Generic wrapper around a user specified combination function, applying that function on a 
 77  * per pixel basis to an image pair, casting the result to an integer, and building an image 
 78  * of the results.  Any additional data required by the combination function can be passed 
 79  * through this function using the void pointer data.  Note that the image resulting from the
 80  * combination is allocated here, and the original image pair are not affected.
 81  */
 82 Imrect         *im_combine(Imrect * im1, Imrect * im2, int *(*func) (), void *data)
 83 {
 84         Imrect         *im3;
 85         Imregion       *roi;
 86         int            *row1, *row2, *row3;
 87         int             lx, ux, ly, uy;
 88         int             width, height;
 89         int             i, j;
 90 
 91         if (im1 == NULL || im2 == NULL)
 92                 return (NULL);
 93 
 94         roi = roi_inter(im1->region, im2->region);
 95         if (roi == NULL)
 96                 return (NULL);
 97         lx = roi->lx;
 98         ux = roi->ux;
 99         ly = roi->ly;
100         uy = roi->uy;
101 
102         width = MIN(im1->width, im2->width);
103         height = MIN(im1->height, im2->height);
104 
105         im3 = im_alloc(height, width, roi, ptr_v);
106         row1 = tvector_alloc(lx, ux, int);
107         row2 = tvector_alloc(lx, ux, int);
108         row3 = tvector_alloc(lx, ux, int);
109 
110         for (i = ly; i < uy; ++i)
111         {
112                 im_get_row(row1, im1, i, lx, ux);
113                 im_get_row(row2, im2, i, lx, ux);
114                 for (j = lx; j < ux; ++j)
115                         row3[j] = (int) (*func) (row1[j], row2[j], data);
116                 im_put_row(row3, im3, i, lx, ux);
117         }
118 
119         tvector_free(row1, lx, void *);
120         tvector_free(row2, lx, void *);
121         tvector_free(row3, lx, void *);
122         rfree(roi);
123         return (im3);
124 }
125 
126 
127 /**
128  * @brief Apply a user-specified per-pixel function to a pair of images to produce a float image.
129  * @param im1 Image to be combined.
130  * @param im2 Image to be combined.
131  * @param func Per-pixel combination function to use.
132  * @param data Pointer to any additional data required by the combination function
133  * @return im3 Integer result of the image combination
134  *
135  * Generic wrapper around a user specified combination function, applying that function on a 
136  * per pixel basis to an image pair, casting the result to a float, and building an image 
137  * of the results.  Any additional data required by the combination function can be passed 
138  * through this function using the void pointer data.  Note that the image resulting from the
139  * combination is allocated here, and the original image pair are not affected.
140  */
141 Imrect         *imf_combine(Imrect * im1, Imrect * im2, double (*func) ( /* ??? */ ), void *data)
142 {
143         Imrect         *im3;
144         Imregion       *roi;
145         float          *row1, *row2, *row3;
146         int             lx, ux, ly, uy;
147         int             width, height;
148         int             i, j;
149 
150         if (im1 == NULL || im2 == NULL)
151                 return (NULL);
152 
153         roi = roi_inter(im1->region, im2->region);
154         if (roi == NULL)
155                 return (NULL);
156         lx = roi->lx;
157         ux = roi->ux;
158         ly = roi->ly;
159         uy = roi->uy;
160 
161         width = MIN(im1->width, im2->width);
162         height = MIN(im1->height, im2->height);
163 
164         im3 = im_alloc(height, width, roi, float_v);
165         row1 = fvector_alloc(lx, ux);
166         row2 = fvector_alloc(lx, ux);
167         row3 = fvector_alloc(lx, ux);
168 
169         for (i = ly; i < uy; ++i)
170         {
171                 im_get_rowf(row1, im1, i, lx, ux);
172                 im_get_rowf(row2, im2, i, lx, ux);
173                 for (j = lx; j < ux; ++j)
174                         row3[j] = (float) (*func) (row1[j], row2[j], data);
175                 im_put_rowf(row3, im3, i, lx, ux);
176         }
177 
178         fvector_free(row1, lx);
179         fvector_free(row2, lx);
180         fvector_free(row3, lx);
181         rfree(roi);
182         return (im3);
183 }
184 
185 
186 
187 /**
188  * @brief Apply a user-specified per-pixel function to a pair of images, writing teh (float) result to im1.
189  * @param im1 Image to be combined.
190  * @param im2 Image to be combined.
191  * @param func Per-pixel combination function to use.
192  * @param data Pointer to any additional data required by the combination function
193  * @return void
194  *
195  * Generic wrapper around a user specified combination function, applying that function on a 
196  * per pixel basis to an image pair, casting the result to an float, and writing the results 
197  * into im1.  Any additional data required by the combination function can be passed 
198  * through this function using the void pointer data.  
199  */
200 void            imf_combine_inplace(Imrect * im1, Imrect * im2, double (*func) ( /* ??? */ ), void *data)
201 {
202         Imregion       *roi;
203         float          *row1, *row2, *row3;
204         int             lx, ux, ly, uy;
205         int             i, j;
206 
207         if (im1 == NULL || im2 == NULL)
208                 return;
209 
210         roi = roi_inter(im1->region, im2->region);
211         rfree(im1->region);
212         im1->region = roi;
213 
214         if (roi == NULL)
215                 return;
216         lx = roi->lx;
217         ux = roi->ux;
218         ly = roi->ly;
219         uy = roi->uy;
220 
221         row1 = fvector_alloc(lx, ux);
222         row2 = fvector_alloc(lx, ux);
223         row3 = fvector_alloc(lx, ux);
224 
225         for (i = ly; i < uy; ++i)
226         {
227                 im_get_rowf(row1, im1, i, lx, ux);
228                 im_get_rowf(row2, im2, i, lx, ux);
229                 for (j = lx; j < ux; ++j)
230                         row3[j] = (float) (*func) (row1[j], row2[j], data);
231                 im_put_rowf(row3, im1, i, lx, ux);
232         }
233 
234         fvector_free(row1, lx);
235         fvector_free(row2, lx);
236         fvector_free(row3, lx);
237         rfree(roi);
238 }
239 
240 
241 /**
242  * @brief Integer pixel-by-pixel sum of an image pair
243  * @param im1 Image to be summed
244  * @param im2 Image to be summed
245  * @return im3 Integer result of the sum.
246  *
247  * Extracts values as integers from an image pair, performs a pixel-by-pixel sum, and writes the result 
248  * to a newly allocated integer image. The sum is performed over the intersection of the regions of the 
249  * original images. This function is intended to be called through im_sum, which selects the relevant 
250  * summing function on the basis of the image types.
251  */
252 Imrect         *imi_sum(Imrect * im1, Imrect * im2)
253 {
254         Imrect         *im3;
255         Imregion       *roi;
256         int            *row1, *row2, *row3;
257         int             lx, ux, ly, uy;
258         int             width, height;
259         int             i, j;
260 
261         if (im1 == NULL || im2 == NULL)
262                 return (NULL);
263 
264         roi = roi_inter(im1->region, im2->region);
265         if (roi == NULL)
266                 return (NULL);
267         lx = roi->lx;
268         ux = roi->ux;
269         ly = roi->ly;
270         uy = roi->uy;
271 
272         width = MIN(im1->width, im2->width);
273         height = MIN(im1->height, im2->height);
274 
275         im3 = im_alloc(height, width, roi, int_v);
276         row1 = ivector_alloc(lx, ux);
277         row2 = ivector_alloc(lx, ux);
278         row3 = ivector_alloc(lx, ux);
279 
280         for (i = ly; i < uy; ++i)
281         {
282                 im_get_row(row1, im1, i, lx, ux);
283                 im_get_row(row2, im2, i, lx, ux);
284                 for (j = lx; j < ux; ++j)
285                         row3[j] = row1[j] + row2[j];
286                 im_put_row(row3, im3, i, lx, ux);
287         }
288 
289         ivector_free(row1, lx);
290         ivector_free(row2, lx);
291         ivector_free(row3, lx);
292         rfree(roi);
293         return (im3);
294 }
295 
296 
297 
298 /**
299  * @brief Float pixel-by-pixel sum of an image pair
300  * @param im1 Image to be summed
301  * @param im2 Image to be summed
302  * @return im3 Float result of the sum.
303  *
304  * Extracts values as floats from an image pair, performs a pixel-by-pixel sum, and writes the result 
305  * to a newly allocated float image. The sum is performed over the intersection of the regions of the 
306  * original images. This function is intended to be called through im_sum, which selects the relevant 
307  * summing function on the basis of the image types.
308  */
309 Imrect         *imf_sum(Imrect * im1, Imrect * im2)
310 {
311         Imrect         *im3;
312         Imregion       *roi;
313         float          *row1, *row2, *row3;
314         int             lx, ux, ly, uy;
315         int             width, height;
316         int             i, j;
317 
318         if (im1 == NULL || im2 == NULL)
319                 return (NULL);
320 
321         roi = roi_inter(im1->region, im2->region);
322         if (roi == NULL)
323                 return (NULL);
324         lx = roi->lx;
325         ux = roi->ux;
326         ly = roi->ly;
327         uy = roi->uy;
328 
329         width = MIN(im1->width, im2->width);
330         height = MIN(im1->height, im2->height);
331 
332         im3 = im_alloc(height, width, roi, float_v);
333         row1 = fvector_alloc(lx, ux);
334         row2 = fvector_alloc(lx, ux);
335         row3 = fvector_alloc(lx, ux);
336 
337         for (i = ly; i < uy; ++i)
338         {
339                 im_get_rowf(row1, im1, i, lx, ux);
340                 im_get_rowf(row2, im2, i, lx, ux);
341                 for (j = lx; j < ux; ++j)
342                         row3[j] = row1[j] + row2[j];
343                 im_put_rowf(row3, im3, i, lx, ux);
344         }
345 
346         fvector_free(row1, lx);
347         fvector_free(row2, lx);
348         fvector_free(row3, lx);
349         rfree(roi);
350         return (im3);
351 }
352 
353 
354 /**
355  * @brief Complex pixel-by-pixel sum of an image pair
356  * @param im1 Image to be summed
357  * @param im2 Image to be summed
358  * @return im3 Complex result of the sum.
359  *
360  * Extracts values as complex numbers from an image pair, performs a pixel-by-pixel sum, and writes the result 
361  * to a newly allocated complex image. The sum is performed over the intersection of the regions of the 
362  * original images. This function is intended to be called through im_sum, which selects the relevant 
363  * summing function on the basis of the image types.
364  */
365 Imrect         *imz_sum(Imrect * im1, Imrect * im2)
366 {
367         Imrect         *im3;
368         Imregion       *roi;
369         Complex        *row1;
370         Complex        *row2;
371         Complex        *row3;
372         int             lx, ux, ly, uy;
373         int             width, height;
374         int             i, j;
375 
376         if (im1 == NULL || im2 == NULL)
377                 return (NULL);
378 
379         roi = roi_inter(im1->region, im2->region);
380         if (roi == NULL)
381                 return (NULL);
382         lx = roi->lx;
383         ux = roi->ux;
384         ly = roi->ly;
385         uy = roi->uy;
386 
387         width = MIN(im1->width, im2->width);
388         height = MIN(im1->height, im2->height);
389 
390         im3 = im_alloc(height, width, roi, complex_v);
391         row1 = zvector_alloc(lx, ux);
392         row2 = zvector_alloc(lx, ux);
393         row3 = zvector_alloc(lx, ux);
394 
395         for (i = ly; i < uy; ++i)
396         {
397                 im_get_rowz(row1, im1, i, lx, ux);
398                 im_get_rowz(row2, im2, i, lx, ux);
399                 for (j = lx; j < ux; ++j)
400                         row3[j] = cmplx_sum(row1[j], row2[j]);
401                 im_put_rowz(row3, im3, i, lx, ux);
402         }
403 
404         zvector_free(row1, lx);
405         zvector_free(row2, lx);
406         zvector_free(row3, lx);
407         rfree(roi);
408         return (im3);
409 }
410 
411 
412 
413 /**
414  * @brief Variable-type sensitive, pixel-by-pixel sum of an image pair.
415  * @param im1 Image to be summed
416  * @param im2 Image to be summed
417  * @return im3 Result of the sum.
418  *
419  * This function wraps imi_sum, imf_sum and imz_sum: it selects the correct summing function on the basis 
420  * of the variable types of the input images and uses it to produce a pixel-by-pixel sum, which is written 
421  * to a newly allocated image and returned.
422  */
423 Imrect         *im_sum(Imrect * im1, Imrect * im2)
424 {
425         if (im1 == NULL || im2 == NULL)
426                 return (NULL);
427         switch (im_sup_vtype(im1->vtype, im2->vtype))
428         {
429         case uchar_v:
430         case char_v:
431         case ushort_v:
432         case short_v:
433         case int_v:
434                 return (imi_sum(im1, im2));
435         case float_v:
436                 return (imf_sum(im1, im2));
437         case complex_v:
438                 return (imz_sum(im1, im2));
439         default:
440                 return (NULL);
441         }
442 }
443 
444 
445 /**
446  * @brief Integer pixel-by-pixel subtraction of an image pair
447  * @param im1 Image to be subtracted from.
448  * @param im2 Image to be subtracted.
449  * @return im3 Integer result of the subtraction.
450  *
451  * Extracts values as integers from an image pair, performs a pixel-by-pixel subtraction, and writes the result 
452  * to a newly allocated integer image. The subtraction is performed over the intersection of the regions of the 
453  * original images. This function is intended to be called through im_diff, which selects the relevant 
454  * subtraction function on the basis of the image types.
455  */
456 Imrect         *imi_diff(Imrect * im1, Imrect * im2)
457 {
458         Imrect         *im3;
459         Imregion       *roi;
460         int            *row1, *row2, *row3;
461         int             lx, ux, ly, uy;
462         int             width, height;
463         int             i, j;
464 
465         if (im1 == NULL || im2 == NULL)
466                 return (NULL);
467 
468         roi = roi_inter(im1->region, im2->region);
469         if (roi == NULL)
470                 return (NULL);
471         lx = roi->lx;
472         ux = roi->ux;
473         ly = roi->ly;
474         uy = roi->uy;
475 
476         width = MIN(im1->width, im2->width);
477         height = MIN(im1->height, im2->height);
478 
479         im3 = im_alloc(height, width, roi, int_v);
480         row1 = ivector_alloc(lx, ux);
481         row2 = ivector_alloc(lx, ux);
482         row3 = ivector_alloc(lx, ux);
483 
484         for (i = ly; i < uy; ++i)
485         {
486                 im_get_row(row1, im1, i, lx, ux);
487                 im_get_row(row2, im2, i, lx, ux);
488                 for (j = lx; j < ux; ++j)
489                         row3[j] = row1[j] - row2[j];
490                 im_put_row(row3, im3, i, lx, ux);
491         }
492 
493         ivector_free(row1, lx);
494         ivector_free(row2, lx);
495         ivector_free(row3, lx);
496         rfree(roi);
497         return (im3);
498 }
499 
500 
501 /**
502  * @brief Float pixel-by-pixel subtraction of an image pair
503  * @param im1 Image to be subtracted from.
504  * @param im2 Image to be subtracted.
505  * @return im3 Float result of the subtraction.
506  *
507  * Extracts values as floats from an image pair, performs a pixel-by-pixel subtraction, and writes the result 
508  * to a newly allocated float image. The subtraction is performed over the intersection of the regions of the 
509  * original images. This function is intended to be called through im_diff, which selects the relevant 
510  * subtraction function on the basis of the image types.
511  */
512 Imrect         *imf_diff(Imrect * im1, Imrect * im2)
513 {
514         Imrect         *im3;
515         Imregion       *roi;
516         float          *row1, *row2, *row3;
517         int             lx, ux, ly, uy;
518         int             width, height;
519         int             i, j;
520 
521         if (im1 == NULL || im2 == NULL)
522                 return (NULL);
523 
524         roi = roi_inter(im1->region, im2->region);
525         if (roi == NULL)
526                 return (NULL);
527         lx = roi->lx;
528         ux = roi->ux;
529         ly = roi->ly;
530         uy = roi->uy;
531 
532         width = MIN(im1->width, im2->width);
533         height = MIN(im1->height, im2->height);
534 
535         im3 = im_alloc(height, width, roi, float_v);
536         row1 = fvector_alloc(lx, ux);
537         row2 = fvector_alloc(lx, ux);
538         row3 = fvector_alloc(lx, ux);
539 
540         for (i = ly; i < uy; ++i)
541         {
542                 im_get_rowf(row1, im1, i, lx, ux);
543                 im_get_rowf(row2, im2, i, lx, ux);
544                 for (j = lx; j < ux; ++j)
545                         row3[j] = row1[j] - row2[j];
546                 im_put_rowf(row3, im3, i, lx, ux);
547         }
548 
549         fvector_free(row1, lx);
550         fvector_free(row2, lx);
551         fvector_free(row3, lx);
552         rfree(roi);
553         return (im3);
554 }
555 
556 
557 /**
558  * @brief Complex pixel-by-pixel subtraction of an image pair
559  * @param im1 Image to be subtracted from.
560  * @param im2 Image to be subtracted.
561  * @return im3 Complex result of the sum.
562  *
563  * Extracts values as complex numbers from an image pair, performs a pixel-by-pixel subtraction, and writes the result 
564  * to a newly allocated complex image. The subtraction is performed over the intersection of the regions of the 
565  * original images. This function is intended to be called through im_diff, which selects the relevant 
566  * subtraction function on the basis of the image types.
567  */
568 Imrect         *imz_diff(Imrect * im1, Imrect * im2)
569 {
570         Imrect         *im3;
571         Imregion       *roi;
572         Complex        *row1;
573         Complex        *row2;
574         Complex        *row3;
575         int             lx, ux, ly, uy;
576         int             width, height;
577         int             i, j;
578 
579         if (im1 == NULL || im2 == NULL)
580                 return (NULL);
581 
582         roi = roi_inter(im1->region, im2->region);
583         if (roi == NULL)
584                 return (NULL);
585         lx = roi->lx;
586         ux = roi->ux;
587         ly = roi->ly;
588         uy = roi->uy;
589 
590         width = MIN(im1->width, im2->width);
591         height = MIN(im1->height, im2->height);
592 
593         im3 = im_alloc(height, width, roi, complex_v);
594         row1 = zvector_alloc(lx, ux);
595         row2 = zvector_alloc(lx, ux);
596         row3 = zvector_alloc(lx, ux);
597 
598         for (i = ly; i < uy; ++i)
599         {
600                 im_get_rowz(row1, im1, i, lx, ux);
601                 im_get_rowz(row2, im2, i, lx, ux);
602                 for (j = lx; j < ux; ++j)
603                         row3[j] = cmplx_diff(row1[j], row2[j]);
604                 im_put_rowz(row3, im3, i, lx, ux);
605         }
606 
607         zvector_free(row1, lx);
608         zvector_free(row2, lx);
609         zvector_free(row3, lx);
610         rfree(roi);
611         return (im3);
612 }
613 
614 
615 /**
616  * @brief Variable-type sensitive, pixel-by-pixel subtraction of an image pair.
617  * @param im1 Image to be subtracted from.
618  * @param im2 Image to be subtracted.
619  * @return im3 Result of the sum.
620  *
621  * This function wraps imi_sum, imf_sum and imz_sum: it selects the correct summing function on the basis 
622  * of the variable types of the input images and uses it to produce a pixel-by-pixel sum, over the
623  * intersection of the input image regions, which is written  to a newly allocated image and returned.
624  */
625 Imrect         *im_diff(Imrect * im1, Imrect * im2)
626 {
627         if (im1 == NULL || im2 == NULL)
628                 return (NULL);
629         switch (im_sup_vtype(im1->vtype, im2->vtype))
630         {
631         case uchar_v:
632         case char_v:
633         case short_v:
634         case ushort_v:
635         case int_v:
636                 return (imi_diff(im1, im2));
637         case float_v:
638                 return (imf_diff(im1, im2));
639         case complex_v:
640                 return (imz_diff(im1, im2));
641         default:
642                 return (NULL);
643         }
644 }
645 
646 
647 /**
648  * @brief Weighted sum of an image pair as a float.
649  * @param a Weighting factor for image 1.
650  * @param b Weighting factor for image 2.
651  * @param im1 Image to be summed.
652  * @param im2 Image to be summed.
653  * @return im3 Float result of teh weighted sum 
654  *
655  * Performs the sum of an image pair, weighted by the factors a and b, on a 
656  * pixel-by-pixel basis over the intersection of the regions of an image pair, 
657  * returning the result as a float image.
658  */
659 Imrect         *imf_wsum(double a, double b, Imrect * im1, Imrect * im2)
660 {
661         Imrect         *im3;
662         Imregion       *roi;
663         float          *row1, *row2, *row3;
664         int             lx, ux, ly, uy;
665         int             width, height;
666         int             i, j;
667 
668         if (im1 == NULL || im2 == NULL)
669                 return (NULL);
670 
671         roi = roi_inter(im1->region, im2->region);
672         if (roi == NULL)
673                 return (NULL);
674         lx = roi->lx;
675         ux = roi->ux;
676         ly = roi->ly;
677         uy = roi->uy;
678 
679         width = MIN(im1->width, im2->width);
680         height = MIN(im1->height, im2->height);
681 
682         im3 = im_alloc(height, width, roi, float_v);
683         row1 = fvector_alloc(lx, ux);
684         row2 = fvector_alloc(lx, ux);
685         row3 = fvector_alloc(lx, ux);
686 
687         for (i = ly; i < uy; ++i)
688         {
689                 im_get_rowf(row1, im1, i, lx, ux);
690                 im_get_rowf(row2, im2, i, lx, ux);
691                 for (j = lx; j < ux; ++j)
692                         row3[j] = (float) (a * row1[j] + b * row2[j]);
693                 im_put_rowf(row3, im3, i, lx, ux);
694         }
695 
696         fvector_free(row1, lx);
697         fvector_free(row2, lx);
698         fvector_free(row3, lx);
699         rfree(roi);
700         return (im3);
701 }
702 
703 
704 /**
705  * @brief Pixel-by-pixel sum of the squares of an image pair.
706  * @param Image to be summed.
707  * @param Image to be summed.
708  * @return Float sum-of-squares image.
709  *
710  * Performs the pixel-by-pixel sum of the squares of an image pair, 
711  * over the intersection of their regions,
712  * writing the result to a newly allocated float image for return.
713  */
714 Imrect         *imf_sumsq(Imrect * im1, Imrect * im2)
715 {
716         Imrect         *im3;
717         Imregion       *roi;
718         float          *row1, *row2, *row3;
719         int             lx, ux, ly, uy;
720         int             width, height;
721         int             i, j;
722 
723         if (im1 == NULL || im2 == NULL)
724                 return (NULL);
725 
726         roi = roi_inter(im1->region, im2->region);
727         if (roi == NULL)
728                 return (NULL);
729         lx = roi->lx;
730         ux = roi->ux;
731         ly = roi->ly;
732         uy = roi->uy;
733 
734         width = MIN(im1->width, im2->width);
735         height = MIN(im1->height, im2->height);
736 
737         im3 = im_alloc(height, width, roi, float_v);
738         row1 = fvector_alloc(lx, ux);
739         row2 = fvector_alloc(lx, ux);
740         row3 = fvector_alloc(lx, ux);
741 
742         for (i = ly; i < uy; ++i)
743         {
744                 im_get_rowf(row1, im1, i, lx, ux);
745                 im_get_rowf(row2, im2, i, lx, ux);
746                 for (j = lx; j < ux; ++j)
747                         row3[j] = row1[j] * row1[j] + row2[j] * row2[j];
748                 im_put_rowf(row3, im3, i, lx, ux);
749         }
750 
751         fvector_free(row1, lx);
752         fvector_free(row2, lx);
753         fvector_free(row3, lx);
754         rfree(roi);
755         return (im3);
756 }
757 
758 
759 /**
760  * @brief Produce an integer image of the maximum pixels, on a pixel-by-pixel basis, from an image pair.
761  * @param im1 Image to be processed.
762  * @param im2 Image to be processed.
763  * @return im3 Integer image of the maximum pixel from each pixel pair.
764  *
765  * Takes an image pair and loops over the pixel positions in the intersection of their regions.  At 
766  * each pixel position, it takes the pixel values from the image pair as integers, selects the largest, 
767  * and writes this to a newly allocated integer image for return. This 
768  * function is intended to be called through im_maxsel, which selects the correct maxsel function on 
769  * the basis of the input image types.
770  */
771 Imrect         *imi_maxsel(Imrect * im1, Imrect * im2)
772 {
773         Imrect         *im3;
774         Imregion       *roi;
775         int            *row1, *row2, *row3;
776         int             lx, ux, ly, uy;
777         int             width, height;
778         int             i, j;
779 
780         if (im1 == NULL || im2 == NULL)
781                 return (NULL);
782 
783         roi = roi_inter(im1->region, im2->region);
784         if (roi == NULL)
785                 return (NULL);
786         lx = roi->lx;
787         ux = roi->ux;
788         ly = roi->ly;
789         uy = roi->uy;
790 
791         width = MIN(im1->width, im2->width);
792         height = MIN(im1->height, im2->height);
793 
794         im3 = im_alloc(height, width, roi, int_v);
795         row1 = ivector_alloc(lx, ux);
796         row2 = ivector_alloc(lx, ux);
797         row3 = ivector_alloc(lx, ux);
798 
799         for (i = ly; i < uy; ++i)
800         {
801                 im_get_row(row1, im1, i, lx, ux);
802                 im_get_row(row2, im2, i, lx, ux);
803                 for (j = lx; j < ux; ++j)
804                         row3[j] = MAX(row1[j], row2[j]);
805                 im_put_row(row3, im3, i, lx, ux);
806         }
807 
808         ivector_free(row1, lx);
809         ivector_free(row2, lx);
810         ivector_free(row3, lx);
811         rfree(roi);
812         return (im3);
813 }
814 
815 
816 
817 /**
818  * @brief Produce a float image of the maximum pixels, on a pixel-by-pixel basis, from an image pair.
819  * @param im1 Image to be processed.
820  * @param im2 Image to be processed.
821  * @return im3 Float image of the maximum pixel from each pixel pair.
822  *
823  * Takes an image pair and loops over the pixel positions in the intersection of their regions.  At 
824  * each pixel position, it takes the pixel values from the image pair as floats, selects the largest, 
825  * and writes this to a newly allocated float image for return. This 
826  * function is intended to be called through im_maxsel, which selects the correct maxsel function on 
827  * the basis of the input image types.
828  */
829 Imrect         *imf_maxsel(Imrect * im1, Imrect * im2)
830 {
831         Imrect         *im3;
832         Imregion       *roi;
833         float          *row1, *row2, *row3;
834         int             lx, ux, ly, uy;
835         int             width, height;
836         int             i, j;
837 
838         if (im1 == NULL || im2 == NULL)
839                 return (NULL);
840 
841         roi = roi_inter(im1->region, im2->region);
842         if (roi == NULL)
843                 return (NULL);
844         lx = roi->lx;
845         ux = roi->ux;
846         ly = roi->ly;
847         uy = roi->uy;
848 
849         width = MIN(im1->width, im2->width);
850         height = MIN(im1->height, im2->height);
851 
852         im3 = im_alloc(height, width, roi, float_v);
853         row1 = fvector_alloc(lx, ux);
854         row2 = fvector_alloc(lx, ux);
855         row3 = fvector_alloc(lx, ux);
856 
857         for (i = ly; i < uy; ++i)
858         {
859                 im_get_rowf(row1, im1, i, lx, ux);
860                 im_get_rowf(row2, im2, i, lx, ux);
861                 for (j = lx; j < ux; ++j)
862                         row3[j] = MAX(row1[j], row2[j]);
863                 im_put_rowf(row3, im3, i, lx, ux);
864         }
865 
866         fvector_free(row1, lx);
867         fvector_free(row2, lx);
868         fvector_free(row3, lx);
869         rfree(roi);
870         return (im3);
871 }
872 
873 
874 /**
875  * @brief Produce a complex image of the maximum pixels, on a pixel-by-pixel basis, from an image pair.
876  * @param im1 Image to be processed.
877  * @param im2 Image to be processed.
878  * @return im3 Complex image of the maximum pixel from each pixel pair.
879  *
880  * Takes an image pair and loops over the pixel positions in the intersection of their regions.  At 
881  * each pixel position, it takes the pixel values from the image pair as complex numbers, selects the 
882  * one with the largest modulus, and writes this to a newly allocated complex image for return. This 
883  * function is intended to be called through im_maxsel, which selects the correct maxsel function on 
884  * the basis of the input image types.
885  */
886 Imrect         *imz_maxsel(Imrect * im1, Imrect * im2)
887 {
888         Imrect         *im3;
889         Imregion       *roi;
890         Complex        *row1;
891         Complex        *row2;
892         Complex        *row3;
893         int             lx, ux, ly, uy;
894         int             width, height;
895         int             i, j;
896 
897         if (im1 == NULL || im2 == NULL)
898                 return (NULL);
899 
900         roi = roi_inter(im1->region, im2->region);
901         if (roi == NULL)
902                 return (NULL);
903         lx = roi->lx;
904         ux = roi->ux;
905         ly = roi->ly;
906         uy = roi->uy;
907 
908         width = MIN(im1->width, im2->width);
909         height = MIN(im1->height, im2->height);
910 
911         im3 = im_alloc(height, width, roi, complex_v);
912         row1 = zvector_alloc(lx, ux);
913         row2 = zvector_alloc(lx, ux);
914         row3 = zvector_alloc(lx, ux);
915 
916         for (i = ly; i < uy; ++i)
917         {
918                 im_get_rowz(row1, im1, i, lx, ux);
919                 im_get_rowz(row2, im2, i, lx, ux);
920                 for (j = lx; j < ux; ++j)
921                 {
922                         if (cmplx_sqrmod(row1[j]) > cmplx_sqrmod(row2[j]))
923                                 row3[j] = row1[j];
924                         else
925                                 row3[j] = row2[j];
926                 }
927                 im_put_rowz(row3, im3, i, lx, ux);
928         }
929 
930         zvector_free(row1, lx);
931         zvector_free(row2, lx);
932         zvector_free(row3, lx);
933         rfree(roi);
934         return (im3);
935 }
936 
937 
938 /**
939  * @brief Produce an image of the maximum pixels, on a pixel-by-pixel basis, from an image pair, with variable type selection.
940  * @param im1 Image to be processed.
941  * @param im2 Image to be processed.
942  * @return im3 Image of the maximum pixel from each pixel pair.
943  *
944  * This function wraps imi_maxsel, imf_maxsel and imz_maxsel, selecting the correct function on the basis 
945  * of the input image variable types. It takes an image pair and loops over the pixel positions in the 
946  * intersection of their regions.  At each pixel position, it takes the pixel values from the image pair,
947  * selects the largest modulus, and writes this to a newly allocated image for return.
948  */
949 Imrect         *im_maxsel(Imrect * im1, Imrect * im2)
950 {
951         if (im1 == NULL || im2 == NULL)
952                 return (NULL);
953         switch (im_sup_vtype(im1->vtype, im2->vtype))
954         {
955         case uchar_v:
956         case char_v:
957         case short_v:
958         case ushort_v:
959         case int_v:
960                 return (imi_maxsel(im1, im2));
961         case float_v:
962                 return (imf_maxsel(im1, im2));
963         case complex_v:
964                 return (imz_maxsel(im1, im2));
965         default:
966                 return (NULL);
967         }
968 }
969 
970 
971 /**
972  * @brief Integer pixel-by-pixel product of an  image pair.
973  * @param im1 Image to be multiplied.
974  * @param im2 Image to be multiplied.
975  * @return im3 Integer image of the pixel-by-pixel product of the input images. 
976  *
977  * Extracts values as integers from an image pair, performs a pixel-by-pixel product, and writes the result 
978  * to a newly allocated integer image. The multiplication is performed over the intersection of the regions of the 
979  * original images. This function is intended to be called through im_prod, which selects the relevant 
980  * multiplication function on the basis of the image types.
981  */
982 Imrect         *imi_prod(Imrect * im1, Imrect * im2)
983 {
984         Imrect         *im3;
985         Imregion       *roi;
986         int            *row1, *row2, *row3;
987         int             lx, ux, ly, uy;
988         int             width, height;
989         int             i, j;
990 
991         if (im1 == NULL || im2 == NULL)
992                 return (NULL);
993 
994         roi = roi_inter(im1->region, im2->region);
995         if (roi == NULL)
996                 return (NULL);
997         lx = roi->lx;
998         ux = roi->ux;
999         ly = roi->ly;
1000         uy = roi->uy;
1001 
1002         width = MIN(im1->width, im2->width);
1003         height = MIN(im1->height, im2->height);
1004 
1005         im3 = im_alloc(height, width, roi, int_v);
1006         row1 = ivector_alloc(lx, ux);
1007         row2 = ivector_alloc(lx, ux);
1008         row3 = ivector_alloc(lx, ux);
1009 
1010         for (i = ly; i < uy; ++i)
1011         {
1012                 im_get_row(row1, im1, i, lx, ux);
1013                 im_get_row(row2, im2, i, lx, ux);
1014                 for (j = lx; j < ux; ++j)
1015                         row3[j] = row1[j] * row2[j];
1016                 im_put_row(row3, im3, i, lx, ux);
1017         }
1018 
1019         ivector_free(row1, lx);
1020         ivector_free(row2, lx);
1021         ivector_free(row3, lx);
1022         rfree(roi);
1023         return (im3);
1024 }
1025 
1026 
1027 /**
1028  * @brief Float pixel-by-pixel product of an  image pair.
1029  * @param im1 Image to be multiplied.
1030  * @param im2 Image to be multiplied.
1031  * @return im3 Float image of the pixel-by-pixel product of the input images. 
1032  *
1033  * Extracts values as floats from an image pair, performs a pixel-by-pixel product, and writes the result 
1034  * to a newly allocated float image. The multiplication is performed over the intersection of the regions of the 
1035  * original images. This function is intended to be called through im_prod, which selects the relevant 
1036  * multiplication function on the basis of the image types.
1037  */
1038 Imrect         *imf_prod(Imrect * im1, Imrect * im2)
1039 {
1040         Imrect         *im3;
1041         Imregion       *roi;
1042         float          *row1, *row2, *row3;
1043         int             lx, ux, ly, uy;
1044         int             width, height;
1045         int             i, j;
1046 
1047         if (im1 == NULL || im2 == NULL)
1048                 return (NULL);
1049 
1050         roi = roi_inter(im1->region, im2->region);
1051         if (roi == NULL)
1052                 return (NULL);
1053         lx = roi->lx;
1054         ux = roi->ux;
1055         ly = roi->ly;
1056         uy = roi->uy;
1057 
1058         width = MIN(im1->width, im2->width);
1059         height = MIN(im1->height, im2->height);
1060 
1061         im3 = im_alloc(height, width, roi, float_v);
1062         row1 = fvector_alloc(lx, ux);
1063         row2 = fvector_alloc(lx, ux);
1064         row3 = fvector_alloc(lx, ux);
1065 
1066         for (i = ly; i < uy; ++i)
1067         {
1068                 im_get_rowf(row1, im1, i, lx, ux);
1069                 im_get_rowf(row2, im2, i, lx, ux);
1070                 for (j = lx; j < ux; ++j)
1071                         row3[j] = row1[j] * row2[j];
1072                 im_put_rowf(row3, im3, i, lx, ux);
1073         }
1074 
1075         fvector_free(row1, lx);
1076         fvector_free(row2, lx);
1077         fvector_free(row3, lx);
1078         rfree(roi);
1079         return (im3);
1080 }
1081 
1082 
1083 /**
1084  * @brief Pixel-by-pixel product of a complex image pair.
1085  * @param im1 Image to be multiplied.
1086  * @param im2 Image to be multiplied.
1087  * @return im3 Complex image of the pixel-by-pixel product of the input images. 
1088  *
1089  * Extracts values from a complex image pair, performs a pixel-by-pixel product, and writes the result 
1090  * to a newly allocated complex image. The multiplication is performed over the intersection of the regions of the 
1091  * original images. This function is intended to be called through im_prod, which selects the relevant 
1092  * multiplication function on the basis of the image types.
1093  */
1094 Imrect         *imz_prod(Imrect * im1, Imrect * im2)
1095 {
1096         Imrect         *im3;
1097         Imregion       *roi;
1098         Complex        *row1;
1099         Complex        *row2;
1100         Complex        *row3;
1101         int             lx, ux, ly, uy;
1102         int             width, height;
1103         int             i, j;
1104 
1105         if (im1 == NULL || im2 == NULL)
1106                 return (NULL);
1107 
1108         roi = roi_inter(im1->region, im2->region);
1109         if (roi == NULL)
1110                 return (NULL);
1111         lx = roi->lx;
1112         ux = roi->ux;
1113         ly = roi->ly;
1114         uy = roi->uy;
1115 
1116         width = MIN(im1->width, im2->width);
1117         height = MIN(im1->height, im2->height);
1118 
1119         im3 = im_alloc(height, width, roi, complex_v);
1120         row1 = zvector_alloc(lx, ux);
1121         row2 = zvector_alloc(lx, ux);
1122         row3 = zvector_alloc(lx, ux);
1123 
1124         for (i = ly; i < uy; ++i)
1125         {
1126                 im_get_rowz(row1, im1, i, lx, ux);
1127                 im_get_rowz(row2, im2, i, lx, ux);
1128                 for (j = lx; j < ux; ++j)
1129                         row3[j] = cmplx_prod(row1[j], row2[j]);
1130                 im_put_rowz(row3, im3, i, lx, ux);
1131         }
1132 
1133         zvector_free(row1, lx);
1134         zvector_free(row2, lx);
1135         zvector_free(row3, lx);
1136         rfree(roi);
1137         return (im3);
1138 }
1139 
1140 
1141 /**
1142  * @brief Pixel-by-pixel product of an image pair, with variable type selection.
1143  * @param im1 Image to be multiplied.
1144  * @param im2 Image to be multiplied.
1145  * @return im3 Image of the pixel-by-pixel product of the input images. 
1146  *
1147  * This function wraps imi_prod, imf_prod and imz_prod, selecting the correct multiplication function on the 
1148  * basis of the input image variable types. It extracts values as from an image pair, performs a pixel-by-pixel 
1149  * product, and writes the result to a newly allocated image. The multiplication is performed over the intersection 
1150  * of the regions of the original images. 
1151  */
1152 Imrect         *im_prod(Imrect * im1, Imrect * im2)
1153 {
1154         if (im1 == NULL || im2 == NULL)
1155                 return (NULL);
1156         switch (im_sup_vtype(im1->vtype, im2->vtype))
1157         {
1158         case uchar_v:
1159         case char_v:
1160         case short_v:
1161         case ushort_v:
1162         case int_v:
1163                 return (imi_prod(im1, im2));
1164         case float_v:
1165                 return (imf_prod(im1, im2));
1166         case complex_v:
1167                 return (imz_prod(im1, im2));
1168         default:
1169                 return (NULL);
1170         }
1171 }
1172 
1173 
1174 /**
1175  * @brief Pixel-by-pixel modal division of a float image pair.
1176  * @param im1 Image acting as numerator.
1177  * @param im2 Image acting as denominator.
1178  * @param thresh Magnitude threshold on the denominator image pixels below which noise stabiliation will be applied.
1179  * @param val Expected maximum value of the denominator distribution?
1180  * @return im3 Float image of the pixel-by-pixel product of the input images. 
1181  *
1182  * Extracts values as floats from an image pair, performs a pixel-by-pixel modal division, and writes the result 
1183  * to a newly allocated float image. The multiplication is performed over the intersection of the regions of the 
1184  * original images. This function is intended to be called through im_div, which selects the relevant 
1185  * multiplication function on the basis of the image types. Modal division is a maximum likelihood technique for 
1186  * stabilising the division against noise effects for pixels where the value of the denominator is similar to its 
1187  * error: see Tina Memo No. 2000-011. 
1188  */
1189 Imrect         *imf_div(Imrect * im1, Imrect * im2, double thresh, double val)
1190 {
1191         Imrect         *im3;
1192         Imregion       *roi;
1193         float          *row1, *row2, *row3;
1194         int             lx, ux, ly, uy;
1195         int             width, height;
1196         int             i, j;
1197 
1198         if (im1 == NULL || im2 == NULL)
1199                 return (NULL);
1200 
1201         roi = roi_inter(im1->region, im2->region);
1202         if (roi == NULL)
1203                 return (NULL);
1204         lx = roi->lx;
1205         ux = roi->ux;
1206         ly = roi->ly;
1207         uy = roi->uy;
1208 
1209         width = MIN(im1->width, im2->width);
1210         height = MIN(im1->height, im2->height);
1211 
1212         im3 = im_alloc(height, width, roi, float_v);
1213         row1 = fvector_alloc(lx, ux);
1214         row2 = fvector_alloc(lx, ux);
1215         row3 = fvector_alloc(lx, ux);
1216 
1217         for (i = ly; i < uy; ++i)
1218         {
1219                 im_get_rowf(row1, im1, i, lx, ux);
1220                 im_get_rowf(row2, im2, i, lx, ux);
1221                 for (j = lx; j < ux; ++j)
1222                 {
1223                         if (row2[j] * row2[j] <= thresh * thresh)
1224                         {
1225                                 if (row2[j] >= 0.0)
1226                                         row3[j] = (float) (row1[j] * 2.0 / (row2[j] + val));
1227                                 else
1228                                         row3[j] = (float) (row1[j] * 2.0 / (row2[j] - val));
1229                         } else
1230                                 row3[j] = row1[j] / row2[j];
1231                 }
1232                 im_put_rowf(row3, im3, i, lx, ux);
1233 
1234         }
1235 
1236         fvector_free(row1, lx);
1237         fvector_free(row2, lx);
1238         fvector_free(row3, lx);
1239         rfree(roi);
1240         return (im3);
1241 }
1242 
1243 
1244 /**
1245  * @brief Pixel-by-pixel modal division of a complex image pair.
1246  * @param im1 Image acting as numerator.
1247  * @param im2 Image acting as denominator.
1248  * @param thresh Magnitude threshold on the denominator image pixels below which noise stabiliation will be applied.
1249  * @param val Expected maximum value of the denominator distribution?
1250  * @return im3 Complex image of the pixel-by-pixel product of the input images. 
1251  *
1252  * Extracts values as complex numbers from an image pair, performs a pixel-by-pixel modal division, and writes the result 
1253  * to a newly allocated complex image. The multiplication is performed over the intersection of the regions of the 
1254  * original images. This function is intended to be called through im_div, which selects the relevant 
1255  * multiplication function on the basis of the image types. Modal division is a maximum likelihood technique for 
1256  * stabilising the division against noise effects for pixels where the value of the denominator is similar to its 
1257  * error: see Tina Memo No. 2000-011. 
1258  */
1259 Imrect         *imz_div(Imrect * im1, Imrect * im2, double thresh, Complex val)
1260 {
1261         Imrect         *im3;
1262         Imregion       *roi;
1263         Complex        *row1;
1264         Complex        *row2;
1265         Complex        *row3;
1266         int             lx, ux, ly, uy;
1267         int             width, height;
1268         int             i, j;
1269 
1270         if (im1 == NULL || im2 == NULL)
1271                 return (NULL);
1272 
1273         roi = roi_inter(im1->region, im2->region);
1274         if (roi == NULL)
1275                 return (NULL);
1276         lx = roi->lx;
1277         ux = roi->ux;
1278         ly = roi->ly;
1279         uy = roi->uy;
1280 
1281         width = MIN(im1->width, im2->width);
1282         height = MIN(im1->height, im2->height);
1283 
1284         im3 = im_alloc(height, width, roi, complex_v);
1285         row1 = zvector_alloc(lx, ux);
1286         row2 = zvector_alloc(lx, ux);
1287         row3 = zvector_alloc(lx, ux);
1288 
1289         for (i = ly; i < uy; ++i)
1290         {
1291                 im_get_rowz(row1, im1, i, lx, ux);
1292                 im_get_rowz(row2, im2, i, lx, ux);
1293                 for (j = lx; j < ux; ++j)
1294                         if (cmplx_sqrmod(row2[j]) <= thresh * thresh)
1295                                 if (row2[j].x > 0.0)
1296                                         row3[j] = cmplx_div(row1[j],
1297                                                             cmplx_times(0.5, cmplx_sum(val, row2[j])));
1298                                 else
1299                                         row3[j] = cmplx_div(row1[j],
1300                                                             cmplx_times(0.5, cmplx_diff(row2[j], val)));
1301                         else
1302                                 row3[j] = cmplx_div(row1[j], row2[j]);
1303                 im_put_rowz(row3, im3, i, lx, ux);
1304         }
1305 
1306         zvector_free(row1, lx);
1307         zvector_free(row2, lx);
1308         zvector_free(row3, lx);
1309         rfree(roi);
1310         return (im3);
1311 }
1312 
1313 
1314 /**
1315  * @brief Pixel-by-pixel modal division of an image pair.
1316  * @param im1 Image acting as numerator.
1317  * @param im2 Image acting as denominator.
1318  * @param thresh Magnitude threshold on the denominator image pixels below which noise stabiliation will be applied.
1319  * @param val Expected maximum value of the denominator distribution?
1320  * @return im3 Image of the pixel-by-pixel product of the input images. 
1321  *
1322  * This function wraps imf_div and imz_dif, selecting the correct model division function on the basis of the input 
1323  * image types. It extracts values from an image pair, performs a pixel-by-pixel modal division, and writes the 
1324  * result to a newly allocated image. The multiplication is performed over the intersection of the regions of the 
1325  * original images. Modal division is a maximum likelihood technique for stabilising the division against noise 
1326  * effects for pixels where the value of the denominator is similar to its error: see Tina Memo No. 2000-011. 
1327  */
1328 Imrect         *im_div(Imrect * im1, Imrect * im2, double thresh, double val)
1329 {
1330         if (im1 == NULL || im2 == NULL)
1331                 return (NULL);
1332         switch (im_sup_vtype(im1->vtype, im2->vtype))
1333         {
1334         case uchar_v:
1335         case char_v:
1336         case short_v:
1337         case ushort_v:
1338         case int_v:
1339         case float_v:
1340                 return (imf_div(im1, im2, thresh, val));
1341         case complex_v:
1342                 return (imz_div(im1, im2, thresh, cmplx(val, 0.0)));
1343         default:
1344                 return (NULL);
1345         }
1346 }
1347 

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