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 /* - hist[] should be declared as having nbin elements
12 - if initialize != 0, sets hist to all zeros,
13 otherwise accumulates into pre-existing histogram counts
14 - values outside the range hbot..htop are not counted
15 */
16
mri_histogram(MRI_IMAGE * im,float hbot,float htop,int initialize,int nbin,int hist[])17 void mri_histogram( MRI_IMAGE *im , float hbot,float htop ,
18 int initialize , int nbin, int hist[] )
19 {
20 register int ih , npix , ii ;
21 register float sbin ;
22 MRI_IMAGE *lim ;
23
24 ENTRY("mri_histogram") ;
25
26 if( im == NULL || htop <= hbot || nbin < 2 ) EXRETURN ;
27
28 /* can handle shorts and floats; all else -> convert datum */
29
30 switch( im->kind ){
31
32 default:
33 STATUS("convert to floats") ;
34 lim = mri_to_float(im) ; break ;
35
36 case MRI_byte:
37 STATUS("convert to shorts") ;
38 lim = mri_to_short(1.0,im) ; break ;
39
40 case MRI_short:
41 case MRI_float:
42 STATUS("keep original") ;
43 lim = im ; break ;
44 }
45
46 npix = lim->nvox ;
47 sbin = 0.999995f * nbin / (htop-hbot) ;
48
49 if( PRINT_TRACING ){
50 char str[256]; sprintf(str,"hbot=%f htop=%f nbin=%d sbin=%f",hbot,htop,nbin,sbin); STATUS(str);
51 }
52
53 if( initialize ) for( ih=0 ; ih < nbin ; ih++ ) hist[ih] = 0 ;
54
55 switch( lim->kind ){
56
57 default: break ;
58
59 case MRI_short:{
60 register short *shar=mri_data_pointer(lim) , val ;
61
62 STATUS("processing shorts") ;
63 for( ii=0 ; ii < npix ; ii++ ){
64 val = shar[ii] ; if( val < hbot || val > htop ) continue ;
65 ih = (int)(sbin*(val-hbot)) ; hist[ih]++ ;
66 }
67 }
68 break ;
69
70 case MRI_float:{
71 register float *flar=mri_data_pointer(lim) , val ;
72
73 STATUS("processing floats") ;
74 for( ii=0 ; ii < npix ; ii++ ){
75 val = flar[ii] ; if( val < hbot || val > htop ) continue ;
76 ih = (int)(sbin*(val-hbot)) ; hist[ih]++ ;
77 }
78 }
79 break ;
80 }
81
82 if( lim != im ) mri_free(lim) ; /* toss temporary array */
83 EXRETURN ;
84 }
85