1 /*****************************************************************************
2 Major portions of this software are copyrighted by the Medical College
3 of Wisconsin, 1994-2000, and are released under the Gnu General Public
4 License, Version 2. See the file README.Copyright for details.
5 ******************************************************************************/
6
7 #include "mrilib.h"
8
9 /*** 7D SAFE ***/
10
11 #ifndef MAX
12 # define MAX(x,y) (((x)>(y))?(x):(y))
13 #endif
14
15 /*---------------------------------------------------------------------------*/
16
mri_edit_image(float pthr,float power,MRI_IMAGE * imin)17 MRI_IMAGE * mri_edit_image( float pthr , float power , MRI_IMAGE *imin )
18 {
19 int ii , npix , nsum ;
20 float val ;
21 MRI_IMAGE *imqq ;
22 float *flin ;
23
24 ENTRY("mri_edit_image") ;
25
26 imqq = mri_to_float( imin ) ;
27 flin = mri_data_pointer( imqq ) ;
28 npix = imqq->nvox ;
29
30 if( (power==0.0 || power==1.0) && (pthr==0.0) ) RETURN(imqq) ;
31
32 if( pthr > 0.0 && pthr < 1.0 ){
33 register float sum , fa , scl,fmax ;
34 register int nsum ;
35
36 fmax = fabs(mri_max(imqq)) ;
37 val = fabs(mri_min(imqq)) ;
38 fmax = MAX(fmax,val) ;
39 val = pthr * fmax ; /* average pixels > pthr * max */
40 sum = 0.0 ;
41 nsum = 0 ;
42
43 for( ii=0 ; ii < npix ; ii++ ){
44 fa = flin[ii] = fabs(flin[ii]) ;
45 if( fa > val ){ sum += fa ; nsum++ ; }
46 }
47 val = pthr * sum / nsum ; /* set threshold based on this */
48
49 #ifdef HARD_THRESH
50 for( ii=0 ; ii < npix ; ii++ ) if(flin[ii] < val) flin[ii] = 0.0 ;
51 #else
52 scl = fmax / (fmax-val) ;
53 for( ii=0 ; ii < npix ; ii++ ){
54 fa = flin[ii] ;
55 flin[ii] = (fa < val) ? (0.0) : (scl*(fa-val)) ;
56 }
57 #endif
58 } /* end of if(pthr) */
59
60 if( power != 0.0 && power != 1.0 ){
61 for( ii=0 ; ii < npix ; ii++ ) flin[ii] = pow( fabs(flin[ii]) , power ) ;
62 }
63
64 MRI_COPY_AUX(imqq,imin) ;
65 RETURN(imqq) ;
66 }
67