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

Linux Cross Reference
Tina6/tina-tools/tinatool/tlmedical/tlmedCoreg_noise.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    :  
 37  * Date    :  
 38  * Version :  
 39  * CVS Id  :  
 40  *
 41  * Notes :
 42  *
 43  *********
 44 */
 45 
 46 #include "tlmedCoreg_noise.h"
 47 
 48 #if HAVE_CONFIG_H
 49 #   include <config.h>
 50 #endif
 51 
 52 
 53 #include <tina/sys/sysDef.h>
 54 #include <tina/sys/sysPro.h>
 55 #include <tina/image/imgDef.h>
 56 #include <tina/image/imgPro.h>
 57 #include <tina/math/mathDef.h>
 58 #include <tina/math/mathPro.h>
 59 #include <tinatool/tlbase/tlbasePro.h>
 60 #include <tinatool/tlmedical/tlmedCoreg_view.h>
 61 
 62 /* Storage for floating image volume copy and new/old version of imptrs */
 63 
 64 static void ***imptrs_store=NULL;
 65 static List *noise_seq=NULL;
 66 
 67 
 68 Imrect *noise_masker(Imrect *im)
 69 {
 70         int j, k;
 71         Imregion *roi=NULL;
 72         Imrect *im1=NULL;
 73 
 74         /* PAB addition: prepare some masks to ensure that noise is not added to zero-padded regions */
 75 
 76         roi = im->region;
 77         im1 = im_alloc((roi->uy-roi->ly), (roi->ux-roi->lx), roi, float_v);
 78         for(j=roi->ly; j<roi->uy; j++)
 79         {
 80                 for(k=roi->lx; k<roi->ux; k++)
 81                 {
 82                         if(im_get_pixf(im, j, k)==0)
 83                         {
 84                                 im_put_pix(0, im1, j, k);
 85                         }
 86                         else
 87                         {
 88                                 im_put_pix(1, im1, j, k);
 89                         }
 90                 }
 91         }
 92         return im1;
 93 }
 94 
 95 
 96 void coreg_add_noise(float noise_multiplier)
 97 {
 98   /* Add noise to the coreg tool's latched images i.e.
 99      the reference images of the coregistration */
100 
101   static Imrect *im_x=NULL, *im_y=NULL, *im_z=NULL;
102   Imrect *im_noise = NULL, *mask_mult=NULL, *immask=NULL;
103   static Imrect *im_sum_x=NULL, *im_sum_y=NULL, *im_sum_z=NULL;
104   Imregion *roi = NULL;
105   double sigma_x, sigma_y;
106   static double ax = 1.0, ay =1.0;
107   static double sigma_noise = 0.0;
108   static int width = 256, height = 256;
109   static int i = 0;
110 
111   if (i == 0)
112   {
113      im_x = get_xcoreg_im();
114      im_y = get_ycoreg_im();
115      im_z = get_zcoreg_im();
116      if(im_z == NULL) return;
117  
118      if ((roi = im_z->region) == NULL) return;
119 
120      sigma_x = imf_diffx_noise(im_z, roi);
121      sigma_y = imf_diffy_noise(im_z, roi);
122      sigma_noise = (sigma_x + sigma_y)/2;
123      format("\nzcoreg Im noise: %f", sigma_noise);
124   }
125 
126   if(im_sum_z!=NULL) im_free(im_sum_z);
127   im_sum_z = NULL;
128   if(im_sum_y!=NULL) im_free(im_sum_y); 
129   im_sum_y = NULL;
130   if(im_sum_x!= NULL) im_free(im_sum_x);
131   im_sum_x = NULL;
132 
133   im_noise = imf_norm_noise(width, height, (int)ax,(int) ay, 0.0, (sigma_noise*noise_multiplier));
134   immask = noise_masker(im_z);
135   mask_mult = im_prod(im_noise, immask);
136   im_sum_z = im_sum(mask_mult,im_z);
137   im_free(immask);
138   im_free(mask_mult);
139   im_free(im_noise);
140 
141   im_noise = imf_norm_noise(width, height, (int)ax,(int) ay, 0.0, (sigma_noise*noise_multiplier));
142   immask = noise_masker(im_y);
143   mask_mult = im_prod(im_noise, immask);
144   im_sum_y = im_sum(mask_mult,im_y);
145   im_free(immask);
146   im_free(mask_mult);
147   im_free(im_noise);
148 
149   im_noise = imf_norm_noise(width, height, (int)ax,(int) ay, 0.0, (sigma_noise*noise_multiplier));
150   immask = noise_masker(im_x);
151   mask_mult = im_prod(im_noise, immask);
152   im_sum_x = im_sum(mask_mult,im_x);
153   im_free(immask);
154   im_free(mask_mult);
155   im_free(im_noise);
156 
157   set_zcoreg_im(im_sum_z);
158   set_ycoreg_im(im_sum_y);
159   set_xcoreg_im(im_sum_x);
160   i++;
161 }
162 
163 
164 void coreg_pab_noise(float noise_multiplier)
165 {
166         /* Add noise to the floating image volume stored as imptrs in
167         coreg_view. */
168 
169         void set_voxel_vtype(Vartype vtype);
170         void ***imptrs_new=NULL;
171         double xcentre, ycentre, zcentre, sigma_x, sigma_y, sigma_noise;
172         Vec3 dummy;
173         Imrect *im_ref_z=NULL, *im_noise=NULL, *im_trans_z=NULL, *imref=NULL, *im_noise_sum;
174         List *get_start_el();
175         List *frame=NULL;
176         Sequence *seq = NULL;
177         int imptrlx, imptrux,imptrly, imptruy, imptrlz, imptruz, i, j, ax=1, ay=1;
178         Imregion *roi=NULL;
179         float pix;
180 
181         /* Check for a double call */
182 
183         if(imptrs_store!=NULL)
184         {
185                 error("Noise addition is being double-called\n", warning);
186                 imptrs_reset();
187         }
188 
189         /* First get the reference images for their roi's */
190 
191         im_ref_z = get_zcoreg_im();
192         coreg_centre(&xcentre, &ycentre, &zcentre, &dummy);
193         if(im_ref_z==NULL) return;
194 
195         /* Now get the floating image z for noise calculations */
196         /* Noise functionality changed here to reflect changes in Imcalc
197         functions: check these */
198         
199         im_trans_z = seq_slicez(zcentre, im_ref_z->region, coreg_bproj);
200         roi = im_trans_z->region;
201         sigma_x = imf_diffx_noise(im_trans_z, roi);
202         sigma_y = imf_diffy_noise(im_trans_z, roi);
203         sigma_noise = (sigma_x + sigma_y)/2;
204         format("\nzcoreg Im noise: %f", sigma_noise);
205 
206         /* Now prepare a new version of imptrs by looping over the sequence
207         and adding noise before writing the rows into imptrs */
208 
209         frame = get_current_seq_start_el();
210         if((frame==NULL)||(imref = frame->to)==NULL) return;
211         roi = imref->region;
212 
213         seq = seq_get_current();
214         imptrlx = roi->lx;
215         imptrux = roi->ux;
216         imptrly = roi->ly;
217         imptruy = roi->uy;
218         imptrlz = seq->offset;
219         imptruz = get_end_frame(seq)+1;
220 
221         for(frame=get_current_seq_start_el(); frame!=NULL; frame=frame->next)
222         {
223                 imref = (Imrect *)frame->to;
224                 roi = imref->region;
225                 im_noise = imf_norm_noise((roi->ux-roi->lx), (roi->uy-roi->ly), ax, ay, 0.0, (sigma_noise*noise_multiplier)); 
226 
227                 for(j=roi->ly; j<roi->uy; j++)
228                 {
229                         for(i=roi->lx; i<roi->ux; i++)
230                         {
231                                 pix = im_get_pixf(imref, j, i);
232                                 if(pix==0.0) im_put_pixf(0.0, im_noise, j, i);
233                         }
234                 }
235 
236                 im_noise_sum = im_sum(im_noise, imref);
237                 noise_seq = list_addtoend(noise_seq, link_alloc(im_noise_sum, 1005));
238                 im_free(im_noise);
239         }
240 
241         imptrs_new = (void ***)parray_alloc(imptrlz,imptrly,imptruz,imptruy);
242         frame = noise_seq;
243         imref = frame->to;
244         for (i=imptrlz;i<imptruz;i++)
245         {
246                 for (j=imptrly;j<imptruy;j++)
247                 {
248                         imptrs_new[i][j] = (void *)IM_ROW(imref,j);
249                 }
250                 if( frame !=NULL)
251                 {
252                         frame = frame->next;
253                         if (frame!=NULL) imref = frame->to;
254                 }
255         } 
256         im_free(im_trans_z);
257 
258         /* Now replace imptrs with the new version */
259 
260         imptrs_store = imptrs_swap(imptrs_new);
261         seq_voxel_vtype(imref->vtype);  
262 }
263 
264 
265 void imptrs_reset()
266 {
267         /* Reset the imptrs in coreg_view, free the new version, and free the noisy images */
268 
269         void set_voxel_vtype(Vartype vtype);
270         void ***imptrs_new=NULL;
271         int imptrlz, imptrly, imptrlx, imptruz, imptruy, imptrux;
272         Imrect *im=NULL;
273         List *frame=NULL;
274         List *get_start_el();
275 
276         if(imptrs_store==NULL) return;
277 
278         imptrs_new = seq_limits(&imptrlz, &imptruz, &imptrly, &imptruy, &imptrlx, &imptrux);
279         imptrs_new = imptrs_swap(imptrs_store);
280 
281         frame = get_current_seq_start_el();
282         im = frame->to;
283         seq_voxel_vtype(im->vtype);     
284 
285         parray_free(imptrs_new,imptrlz,imptrly,imptruz,imptruy);
286         imptrs_new=NULL;
287 /*
288         for (frame = noise_seq; frame!=NULL; frame = frame->next)
289         {
290                 im = (Imrect *)frame->to;               
291                 im_free(im);
292                 rfree((void *)frame);
293         }
294         noise_seq=NULL;
295 */
296         dd_list_rm(noise_seq, im_free);
297         noise_seq = NULL;
298         imptrs_store=NULL;
299 }
300 

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