Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
spanner.c
Go to the documentation of this file.00001 #include "mrilib.h"
00002
00003 typedef struct {
00004 int bot , top ;
00005 float pval ;
00006 } spanfit ;
00007
00008 float evaluate_span( int , int , int , int , float * , float ** ) ;
00009
00010
00011
00012 spanfit find_best_span( int ndim , int nvec , int minspan ,
00013 float * cvec , float ** bvec )
00014 {
00015 spanfit result = {0,0,0.0} ;
00016 int ii,kk , bot,top , bot_best,top_best ;
00017 float val , val_best ;
00018
00019 if( minspan < 3 || ndim < minspan || nvec < 100 ) return result ;
00020 if( cvec == NULL || bvec == NULL ) return result ;
00021
00022 val_best = -1.0 ;
00023 for( bot=0 ; bot < ndim+1-minspan ; bot++ ){
00024 for( top=bot+minspan-1 ; top < ndim ; top++ ){
00025 val = evaluate_span( ndim,nvec , bot,top , cvec,bvec ) ;
00026 if( val > val_best ){
00027 val_best = val ; bot_best = bot ; top_best = top ;
00028 }
00029 }
00030 }
00031
00032 evaluate_span( 0,0,0,0,NULL,NULL ) ;
00033
00034 result.bot = bot_best ; result.top = top_best ; result.pval = val_best ;
00035 return result ;
00036 }
00037
00038
00039
00040 float evaluate_span( int ndim, int nvec,
00041 int bot , int top , float * cvec , float ** bvec )
00042 {
00043 int kk , ibot=bot,itop=top , nneg ;
00044 float cbar, *qvec ;
00045 register int ii ;
00046 register float sum ;
00047
00048 static float * bsum=NULL , * cnorm=NULL ;
00049 static int nbsum=0 , ncnorm=0 ;
00050
00051 if( nvec > nbsum ){
00052 if( bsum != NULL ) free(bsum) ;
00053 bsum = (float *) malloc(sizeof(float)*nvec) ;
00054 nbsum = nvec ;
00055 } else if( nvec <= 0 ){
00056 if( bsum != NULL ){ free(bsum) ; bsum = NULL; nbsum = 0; }
00057 if( cnorm != NULL ){ free(cnorm); cnorm = NULL; ncnorm = 0; }
00058 return 0.0 ;
00059 }
00060
00061 if( ndim > ncnorm ){
00062 if( cnorm != NULL ) free(cnorm) ;
00063 cnorm = (float *) malloc(sizeof(float)*ndim) ;
00064 ncnorm = ndim ;
00065 }
00066
00067
00068
00069 sum = 0.0 ;
00070 for( ii=ibot ; ii <= itop ; ii++ ) sum += cvec[ii] ;
00071 cbar = sum/(itop-ibot+1) ; sum = 0.0 ;
00072 for( ii=ibot ; ii <= itop ; ii++ ){
00073 cnorm[ii] = cvec[ii] - cbar ; sum += cnorm[ii]*cnorm[ii] ;
00074 }
00075 if( sum <= 0.0 ) return 1.0 ;
00076
00077
00078
00079 for( kk=0 ; kk < nvec ; kk++ ){
00080 qvec = bvec[kk] ; sum = 0.0 ;
00081 for( ii=ibot ; ii <= itop ; ii++ ) sum += qvec[ii] * cnorm[ii] ;
00082 bsum[kk] = sum ;
00083 }
00084
00085
00086
00087 for( nneg=ii=0 ; ii < nvec ; ii++ ) if( bsum[ii] <= 0.0 ) nneg++ ;
00088
00089
00090
00091 sum = (float)nneg / (float)nvec ;
00092 return sum ;
00093 }
