Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
mri_blur3d.c
Go to the documentation of this file.00001 #include "mrilib.h"
00002
00003
00004
00005 #define WT0 0.3
00006 #define WT1 0.2
00007 #define WT2 0.15
00008
00009 void MRI_5blur_inplace_3D( MRI_IMAGE *im )
00010 {
00011 int ii,jj,kk , nx,ny,nz , nxy,nxyz , off,top ;
00012 float *f , *r ;
00013
00014 ENTRY("MRI_5blur_inplace_3D") ;
00015
00016 if( im == NULL || im->kind != MRI_float ) EXRETURN ;
00017
00018 nx = im->nx ; ny = im->ny ; nz = im->nz ; nxy = nx*ny ; nxyz = nxy*nz ;
00019 f = MRI_FLOAT_PTR(im) ;
00020
00021
00022
00023 #undef D
00024 #define D 1
00025 if( nx > 3 ){
00026 top = nx-2 ;
00027 r = (float *) malloc(sizeof(float)*nx) ;
00028 for( kk=0 ; kk < nz ; kk++ ){
00029 for( jj=0 ; jj < ny ; jj++ ){
00030 off = jj*nx + kk*nxy ;
00031 r[0] = WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00032 r[1] = WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00033 for( ii=2 ; ii < top ; ii++,off+=D ){
00034 r[ii] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D];
00035 }
00036 r[ii] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]; off+=D; ii++;
00037 r[ii] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off] ;
00038 off = jj*nx + kk*nxy ;
00039 for( ii=0 ; ii < nx ; ii++ ) f[ii*D+off] = r[ii] ;
00040 }
00041 }
00042 free(r) ;
00043 }
00044
00045
00046
00047 #undef D
00048 #define D nx
00049 if( ny > 3 ){
00050 top = ny-2 ;
00051 r = (float *) malloc(sizeof(float)*ny) ;
00052 for( kk=0 ; kk < nz ; kk++ ){
00053 for( ii=0 ; ii < nx ; ii++ ){
00054 off = ii + kk*nxy ;
00055 r[0] = WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00056 r[1] = WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00057 for( jj=2 ; jj < top ; jj++,off+=D ){
00058 r[jj] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D];
00059 }
00060 r[jj] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]; off+=D; jj++;
00061 r[jj] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off] ;
00062 off = ii + kk*nxy ;
00063 for( jj=0 ; jj < ny ; jj++ ) f[jj*D+off] = r[jj] ;
00064 }
00065 }
00066 free(r) ;
00067 }
00068
00069
00070
00071 #undef D
00072 #define D nxy
00073 if( nz > 3 ){
00074 top = nz-2 ;
00075 r = (float *) malloc(sizeof(float)*nz) ;
00076 for( jj=0 ; jj < ny ; jj++ ){
00077 for( ii=0 ; ii < nx ; ii++ ){
00078 off = ii + jj*nx ;
00079 r[0] = WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00080 r[1] = WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D]; off+=D;
00081 for( kk=2 ; kk < top ; kk++,off+=D ){
00082 r[kk] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]+WT2*f[off+2*D];
00083 }
00084 r[kk] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off]+WT1*f[off+D]; off+=D; kk++;
00085 r[kk] = WT2*f[off-2*D]+WT1*f[off-D]+WT0*f[off] ;
00086 off = ii + jj*nx ;
00087 for( kk=0 ; kk < nz ; kk++ ) f[kk*D+off] = r[kk] ;
00088 }
00089 }
00090 free(r) ;
00091 }
00092
00093 EXRETURN ;
00094 }