Doxygen Source Code Documentation
mri_stat_seq.c File Reference
#include "mrilib.h"
Go to the source code of this file.
Functions | |
MRI_IMAGE ** | mri_stat_seq (MRI_IMAGE *imin) |
Function Documentation
|
Definition at line 18 of file mri_stat_seq.c. References EXIT, MRI_IMAGE::kind, malloc, MRI_COPY_AUX, mri_data_pointer(), mri_free(), MRI_IS_2D, mri_new(), mri_to_float(), MRI_IMAGE::nx, MRI_IMAGE::ny, and retval(). Referenced by main(), and mri_align_crao().
00019 { 00020 00021 /* static variables (exist between calls) */ 00022 00023 static MRI_IMAGE * imsum , * imsumq ; 00024 static int nim = 0 , npix , nx , ny ; 00025 static double * dbsum , * dbsumq ; 00026 00027 /* local variables */ 00028 00029 register int ii ; 00030 MRI_IMAGE * imfl , * imsd , ** retval ; 00031 float * flar , * sdar ; 00032 double scl , vscl , vvv ; 00033 00034 /*** case: set up new problem ***/ 00035 00036 if( nim == 0 ){ 00037 00038 if( imin == NULL ){ 00039 fprintf(stderr,"mri_stat_seq: NULL argument for initial call!\n") ; 00040 EXIT(1) ; 00041 } 00042 00043 if( ! MRI_IS_2D(imin) ){ 00044 fprintf(stderr,"\n*** mri_stat_seq: only works on 2D images!\n") ; 00045 EXIT(1) ; 00046 } 00047 00048 nx = imin->nx ; 00049 ny = imin->ny ; 00050 npix = nx * ny ; 00051 00052 imsum = mri_new( nx , ny , MRI_double ) ; 00053 imsumq = mri_new( nx , ny , MRI_double ) ; 00054 dbsum = mri_data_pointer( imsum ) ; 00055 dbsumq = mri_data_pointer( imsumq ) ; 00056 00057 MRI_COPY_AUX(imsum,imin) ; 00058 00059 for( ii=0 ; ii < npix ; ii++ ) dbsum[ii] = dbsumq[ii] = 0.0 ; 00060 00061 } 00062 00063 /*** case: add new image into problem ***/ 00064 00065 if( imin != NULL ){ 00066 00067 if( imin->nx != nx || imin->ny != ny ){ 00068 fprintf(stderr,"mri_stat_seq: input image size mismatch!\n") ; 00069 EXIT(1) ; 00070 } 00071 00072 if( ! MRI_IS_2D(imin) ){ 00073 fprintf(stderr,"\n*** mri_stat_seq: only works on 2D images!\n") ; 00074 EXIT(1) ; 00075 } 00076 00077 if( imin->kind == MRI_float ){ 00078 imfl = imin ; 00079 } else { 00080 imfl = mri_to_float( imin ) ; 00081 } 00082 flar = mri_data_pointer( imfl ) ; 00083 00084 for( ii=0 ; ii < npix ; ii++ ){ 00085 dbsum[ii] += flar[ii] ; 00086 dbsumq[ii] += flar[ii] * flar[ii] ; 00087 } 00088 00089 if( imin != imfl ) mri_free( imfl ) ; 00090 nim++ ; 00091 return NULL ; 00092 } 00093 00094 /*** case: make report on results, and reset static data ***/ 00095 00096 if( nim < 1 ){ 00097 fprintf(stderr,"mri_stat_seq: # images input=%d; too small!\n",nim) ; 00098 EXIT(1) ; 00099 } 00100 00101 imfl = mri_new( nx , ny , MRI_float ) ; 00102 imsd = mri_new( nx , ny , MRI_float ) ; 00103 flar = mri_data_pointer( imfl ) ; 00104 sdar = mri_data_pointer( imsd ) ; 00105 00106 MRI_COPY_AUX(imfl,imsum) ; 00107 00108 scl = 1.0 / nim ; 00109 vscl = (nim==1) ? 1.0 : sqrt( ((double)(nim))/(nim-1.0) ) ; 00110 00111 for( ii=0 ; ii < npix ; ii++ ){ 00112 flar[ii] = scl * dbsum[ii] ; 00113 vvv = scl * dbsumq[ii] - flar[ii] * flar[ii] ; 00114 sdar[ii] = (vvv > 0.0) ? (sqrt(vvv)*vscl) : 0.0 ; 00115 } 00116 00117 mri_free( imsum ) ; mri_free( imsumq ) ; 00118 nim = 0 ; 00119 00120 retval = (MRI_IMAGE **) malloc( 2 * sizeof(MRI_IMAGE *) ) ; 00121 if( retval == NULL ){ 00122 fprintf(stderr,"mri_stat_seq: malloc error for retval!\n") ; 00123 EXIT(1) ; 00124 } 00125 00126 retval[0] = imfl ; 00127 retval[1] = imsd ; 00128 00129 return retval ; 00130 } |