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

Linux Cross Reference
Tina5/tina-tools/tinatool/tlmedical/tlmedSeg_hist.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-tools/tinatool/tlmedical/tlmedSeg_hist.c,v $
 37  * Date    :  $Date: 2005/02/07 23:29:19 $
 38  * Version :  $Revision: 1.5 $
 39  * CVS Id  :  $Id: tlmedSeg_hist.c,v 1.5 2005/02/07 23:29:19 paul Exp $
 40  *
 41  * Notes :
 42  *
 43  *********
 44 */
 45 
 46 #include "tlmedSeg_hist.h"
 47 
 48 #if HAVE_CONFIG_H
 49 #   include <config.h>
 50 #endif
 51 
 52 #include <stdio.h>
 53 #include <math.h>
 54 #include <tina/sys/sysDef.h>
 55 #include <tina/sys/sysPro.h>
 56 #include <tina/math/mathDef.h>
 57 #include <tina/math/mathPro.h>
 58 #include <tina/image/imgDef.h>
 59 #include <tina/image/imgPro.h>
 60 #include <tinatool/draw/drawDef.h>
 61 #include <tinatool/draw/drawPro.h>
 62 #include <tinatool/tlbase/tlbasePro.h>
 63 
 64 
 65 #define GMAX 7.0
 66 
 67 /*
 68 double pure_density_plot(Mixmodel *m, float *x, int blob, int k)
 69 {   
 70    double z, tmp;
 71 
 72    tmp = (x[k]-m->vectors[blob][k])*m->cov[blob][k][k];   
 73    z = tmp*(x[k]-m->vectors[blob][k]);
 74    if (fabs(z)>GMAX*GMAX) return (0.0);   
 75    z  = exp(-z/2.0)/sqrt(TWOPI/m->cov[blob][k][k]);  
 76    return(z);
 77 }
 78 */
 79 
 80 double em_graph_data(Mixmodel *m, float pix_val)
 81 {
 82    int i,j; 
 83    double tot = 0.0, a =0.0, b;
 84    float *pix_val_arr;
 85  
 86    pix_val_arr = (&pix_val); 
 87 
 88    for( i = 0; i < m->nmix; i++)
 89    {
 90       for(j = 0; j < m->nmix; j++)
 91       { 
 92          if(i == j)
 93             a = m->normal[i]*pure_density(m,pix_val_arr,i);              
 94          else
 95          {
 96             b = part_fraction(m,pix_val_arr,i,j);
 97             a = m->par_dens[i][j]*part_density(m,pix_val_arr,i,j,b);
 98          }
 99          tot += a;
100       }    
101    }
102    tot += m->background;
103    return(tot);  
104 }
105 
106 static double em_graph_data_pure(Mixmodel *m, float pix_val)
107 {
108    int i; 
109    double tot = 0.0, a =0.0;
110    float *pix_val_arr;
111  
112    pix_val_arr = (&pix_val); 
113 
114 
115    for( i = 0; i < m->nmix; i++) 
116    {
117       a = m->normal[i]*pure_density(m,pix_val_arr,i);              
118       tot += a;          
119    }
120    tot += m->background;
121    return(tot);  
122 }
123 
124 
125 static double em_graph_data_mix(Mixmodel *m, float pix_val)
126 {
127    int i,j; 
128    double tot = 0.0, a =0.0,b;
129    float *pix_val_arr;
130 
131    pix_val_arr = (&pix_val); 
132 
133    for( i = 0; i < m->nmix; i++)
134    {
135       for(j = 0; j < m->nmix; j++)
136       { 
137          if(i != j)
138          {
139             b = part_fraction(m,pix_val_arr,i,j);
140             a = m->par_dens[i][j]*part_density(m,pix_val_arr,i,j,b);
141             tot += a;
142          }
143       }    
144    }
145 /*
146    tot += model->background;
147 */
148    return(tot);  
149 }
150 
151 
152 static void em_graph_fit(Tv *tv, shistogram *hist, int im_no, Mixmodel *m)
153 {
154    int i, ndata;
155    float *xdata,*ydata, *zdata,*zdata1,*zdata2;
156 
157    m = get_submodel(m, im_no, im_no);
158    if (tv == NULL || tv->tv_screen == NULL || m==NULL) return;
159    ndata = 2*hist->xbins;
160    xdata = fvector_alloc(0, ndata);
161    ydata = fvector_alloc(0, ndata);
162    zdata = fvector_alloc(0, ndata);
163    zdata1 = fvector_alloc(0, ndata);
164    zdata2 = fvector_alloc(0, ndata);
165   
166    for (i = 0;i < hist->xbins;i++)
167    {
168       xdata[2*i] =  (float)(i+0.0001)*hist->xincr+hist->xmin;
169       xdata[2*i+1] =  (float)(i+1.0001)*hist->xincr+hist->xmin;
170       ydata[2*i] =  hfill1(hist,xdata[2*i],0.0);
171       ydata[2*i+1] =  hfill1(hist,xdata[2*i],0.0);    
172       zdata[2*i] = hist->xincr*(float) em_graph_data(m, xdata[2*i]);
173       zdata[2*i+1] = hist->xincr*(float) em_graph_data(m, xdata[2*i+1]);
174       zdata1[2*i] = hist->xincr*(float) em_graph_data_pure(m, xdata[2*i]);
175       zdata1[2*i+1] = hist->xincr*(float) em_graph_data_pure(m, xdata[2*i+1]);
176       zdata2[2*i] = hist->xincr*(float) em_graph_data_mix(m, xdata[2*i]);
177       zdata2[2*i+1] = hist->xincr*(float) em_graph_data_mix(m, xdata[2*i+1]);
178    }
179 
180    tv_erase(tv);
181 
182    plot(PL_INIT, PL_TV, tv,
183         PL_AXIS_COLOR, black,
184         PL_TITLE, "Histogram Fit",
185         PL_COLOR, red,
186         PL_GRAPH_DATA, (int) ndata, (float*)xdata, (float*)ydata,
187         PL_COLOR, blue,
188         PL_GRAPH_DATA, (int) ndata, (float*)xdata, (float*)zdata,
189         PL_COLOR, green,
190         PL_GRAPH_DATA, (int) ndata, (float*)xdata, (float*)zdata1,
191         PL_COLOR, magenta,
192         PL_GRAPH_DATA, (int) ndata, (float*)xdata, (float*)zdata2,
193         PL_PLOT,
194         NULL);
195 
196    fvector_free(xdata,0);
197    fvector_free(ydata,0);
198    fvector_free(zdata,0);
199    fvector_free(zdata1,0);
200    fvector_free(zdata2,0);
201    mixmodel_free(m);
202 }
203 
204 void  em_hist_plot(Mixmodel *m)
205 {
206     List *get_seq_current_el(Sequence *seq);
207     List *lptr;
208     Imrect *im = NULL, *histim;
209     static shistogram *hist = NULL;
210     Imrect *im_hist(shistogram **imhist, double k, double range,
211                      Imrect *im, Imregion *roi);
212     double range,k; /* range around mean; k - histogram mean */
213                     /* hist->xmin = k-range, hist->xmax= k+range */
214     int  im_no,xmax,ymax,xmin,ymin;
215     char temp[1024];
216     double max_g,min_g;
217 
218     Imregion *roi = tv_get_im_roi(seq_tv_get());   
219     Sequence *seq = seq_get_current();
220    
221     if(seq == NULL)
222     {
223        error("There is no sequence! \n", warning);
224        return; 
225     }
226     if(m == NULL)
227     {
228        error("Model parameters not defined!\n", warning);
229        return;
230     }
231     lptr = get_seq_current_el(seq);
232     im = (Imrect *)(lptr->to); 
233     if (im == NULL)
234     {
235        error("Image sequence is missing!", warning);
236        return;
237     }
238     im_no = get_seq_current_frame(seq);
239 
240     /* Check if number of images in sequence are more than number of parameters specified */
241     if (im_no  >= m->ndim) 
242     {
243        error("More images than parameters - using the first params ", warning);
244        im_no = 0;
245     }
246     max_g = imf_max(im, &ymax,&xmax);
247     min_g = imf_min(im, &ymin,&xmin);
248     k = (min_g + max_g)/2;
249     range = fabs(max_g-min_g)/2; 
250 
251     histim = im_hist(&hist, k, range, im, roi);
252     im_free(histim);
253 /*
254     stack_push(histim, IMRECT, im_free);
255 */
256     sprintf(temp, "overflows %f underflows %f entries %f mean %f min %f max %f\n", hist->over,hist->under,(float)hist->entries,(float)hist->mean, hist->xmin, hist->xmax);
257 /*
258     format(temp);
259 */
260     em_graph_fit((Tv *)imcalc_graph_tv_get(),hist,im_no,m);
261     if (hist != NULL)
262        hfree(hist);
263     hist = NULL;   
264 }
265 

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