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

Linux Cross Reference
Tina6/tina-libs/tina/image/imgPrc_erf.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_erf.c,v $
 37  * Date    :  $Date: 2003/09/30 16:53:49 $
 38  * Version :  $Revision: 1.6 $
 39  * CVS Id  :  $Id: imgPrc_erf.c,v 1.6 2003/09/30 16:53:49 ipoole Exp $
 40  *
 41  * Author  : Legacy TINA
 42  */
 43 
 44 /** 
 45  *  @file
 46  *  @brief Apply the erf() (cummulative normal error) function to image values, including 
 47  *  complex types. 
 48  *
 49  *  Useful for comuting hypothesis probabilities and testing the validity
 50  *  of Gaussian assumtions in data distributions (erf is flat for a Gaussian).
 51  *
 52  */
 53 
 54 #include "imgPrc_erf.h"
 55 
 56 #if HAVE_CONFIG_H
 57 #include <config.h>
 58 #endif
 59 
 60 #include <math.h>
 61 #include <float.h>
 62 #include <limits.h>
 63 #include <tina/sys/sysDef.h>
 64 #include <tina/sys/sysPro.h>
 65 #include <tina/math/mathDef.h>
 66 #include <tina/math/mathPro.h>
 67 #include <tina/image/img_GenDef.h>
 68 #include <tina/image/img_GenPro.h>
 69 
 70 // why is thresh passed by reference - it seems to be an input parameter only?  ipoole
 71 Imrect         *imf_erf(Imrect * im1, double *thresh, double *ml)
 72  
 73 {
 74         Imrect         *im2;
 75         Imregion       *roi;
 76         float          *row1, *row2;
 77         int             lx, ux, ly, uy;
 78         int             i, j;
 79 
 80         if (im1 == NULL)
 81                 return (NULL);
 82 
 83         roi = im1->region;
 84         if (roi == NULL)
 85                 return (NULL);
 86         lx = roi->lx;
 87         ux = roi->ux;
 88         ly = roi->ly;
 89         uy = roi->uy;
 90 
 91         im2 = im_alloc(im1->height, im1->width, roi, float_v);
 92         row1 = fvector_alloc(lx, ux);
 93         row2 = fvector_alloc(lx, ux);
 94 
 95         for (i = ly; i < uy; ++i)
 96         {
 97                 im_get_rowf(row1, im1, i, lx, ux);
 98                 for (j = lx; j < ux; ++j)
 99                 {
100                         row2[j] = (float) 0.50 *(erfc(0.7071 * (*thresh - row1[j])));
101                         if (row2[j] <= sqrt(FLT_MIN))
102                                 row2[j] = sqrt(FLT_MIN);
103                 }
104                 im_put_rowf(row2, im2, i, lx, ux);
105         }
106 
107         fvector_free(row1, lx);
108         fvector_free(row2, lx);
109         *ml = 0.0;
110         return (im2);
111 }
112 
113 Imrect         *imz_erf(Imrect * im1, double *thresh, double *ml)
114 {
115         Imrect         *im2;
116         Imregion       *roi;
117         Complex        *row1;
118         float          *row2;
119         char            temp[260];
120         double          ml_total = 0.0;
121         int             vox_tot = 0.0;
122         int             lx, ux, ly, uy;
123         int             i, j;
124 
125         if (im1 == NULL)
126                 return (NULL);
127 
128         roi = im1->region;
129         if (roi == NULL)
130                 return (NULL);
131         lx = roi->lx;
132         ux = roi->ux;
133         ly = roi->ly;
134         uy = roi->uy;
135 
136         im2 = im_alloc(im1->height, im1->width, roi, float_v);
137         row1 = zvector_alloc(lx, ux);
138         row2 = fvector_alloc(lx, ux);
139 
140         ml_total = 0.0;
141         for (i = ly; i < uy; ++i)
142         {
143                 im_get_rowz(row1, im1, i, lx, ux);
144                 for (j = lx; j < ux; ++j)
145                 {
146                         if (cmplx_re(row1[j]) == 1)
147                         {
148                                 vox_tot++;
149                                 ml_total += cmplx_im(row1[j]);
150                         }
151                 }
152         }
153         sprintf(temp, "total z %f ndf %d \n", ml_total / sqrt(vox_tot), vox_tot);
154         format(temp);
155         for (i = ly; i < uy; ++i)
156         {
157                 im_get_rowz(row1, im1, i, lx, ux);
158                 for (j = lx; j < ux; ++j)
159                 {
160                         if (cmplx_re(row1[j]) == 1)
161                         {
162                                 row2[j] = (float) 0.5 *(erfc(0.7071 * (*thresh) - 0.5 * cmplx_im(row1[j])));
163                                 if (row2[j] <= sqrt(FLT_MIN))
164                                         row2[j] = sqrt(FLT_MIN);
165                         } else if (cmplx_re(row1[j]) == 0)
166                         {
167                                 row2[j] = 1.0 - 0.5 * (erfc(0.7071 * (*thresh) - 0.5 * cmplx_im(row1[j])));
168                                 if (row2[j] <= sqrt(FLT_MIN))
169                                         row2[j] = sqrt(FLT_MIN);
170                         }
171                 }
172                 im_put_rowf(row2, im2, i, lx, ux);
173         }
174 
175 
176         zvector_free(row1, lx);
177         fvector_free(row2, lx);
178         *ml = ml_total;
179         return (im2);
180 }
181 
182 Imrect         *im_erf(Imrect * im, double *thresh, double *ml)
183 {
184         switch (im->vtype)
185         {
186                 case complex_v:
187                 return (imz_erf(im, thresh, ml));
188         default:
189                 return (imf_erf(im, thresh, ml));
190         }
191 }
192 

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