00001 /*****************************************************************************
00002    Major portions of this software are copyrighted by the Medical College
00003    of Wisconsin, 1994-2000, and are released under the Gnu General Public
00004    License, Version 2.  See the file README.Copyright for details.
00005 ******************************************************************************/
00006    
00007 #include "mrilib.h"
00008 
00009 MRI_IMAGE * ts_to_ftime( int nwin , MRI_IMAGE * tim )
00010 {
00011    MRI_IMAGE * ftim ;
00012    float * ftar , * tar , * xar , * win , val ;
00013    complex * cxar ;
00014    int nt , nf , nw2 , nfft, nfft2 , ii,it,kf ;
00015 
00016    if( tim == NULL || tim->kind != MRI_float ||
00017        tim->ny > 1 || nwin < 16              || tim->nx < nwin ) return NULL ;
00018 
00019    nt = tim->nx ; tar = MRI_FLOAT_PTR(tim) ;
00020    if( nwin%2 == 0 ) nwin++ ;
00021    nw2 = nwin/2 ;
00022 
00023    xar = (float *) malloc(sizeof(float)*(nt+2*nw2)) ;
00024    for( ii=0 ; ii < nw2 ; ii++ ) xar[ii] = xar[nt+nw2+ii] = 0.0 ;
00025    memcpy( xar+nw2 , tar , sizeof(float)*nt ) ;
00026 
00027    nfft = csfft_nextup_one35( nwin ) ; nfft2 = nfft/2 ;
00028    ftim = mri_new( nt , nfft2 , MRI_float ) ;
00029    ftar = MRI_FLOAT_PTR(ftim) ;
00030    cxar = (complex *) malloc( sizeof(complex) * nfft ) ;
00031 
00032 fprintf(stderr,"nfft = %d\n",nfft) ;
00033 
00034    win = mri_setup_taper( nwin , 1.0 ) ;
00035 
00036    for( it=0 ; it < nt ; it++ ){
00037 
00038       for( val=0.0,ii=0 ; ii < nwin ; ii++ ) val += xar[it+ii-nw2] ;
00039       val /= nwin ;
00040 
00041       for( ii=0 ; ii < nwin ; ii++ ){
00042          cxar[ii].r = (xar[it+ii-nw2]-val)*win[ii] ; cxar[ii].i = 0.0 ;
00043       }
00044       for( ii=nwin ; ii < nfft ; ii++ ) cxar[ii].r = cxar[ii].i = 0.0 ;
00045 
00046       csfft_cox( -1 , nfft , cxar ) ;
00047 
00048       for( kf=0 ; kf < nfft2 ; kf++ )
00049          ftar[it+nt*kf] = CSQR( cxar[kf] ) ;
00050    }
00051 
00052    free(win) ; free(xar) ; free(cxar) ;
00053    return ftim ;
00054 }
00055 
00056 int main( int argc , char * argv[] )
00057 {
00058    MRI_IMAGE * inim , * outim ;
00059    int nopt,nwin;
00060 
00061    /*-- help? --*/
00062 
00063    if( argc < 4 || strcmp(argv[1],"-help") == 0 ){
00064      printf("Usage: 1dftime nwin infile outfile\n"
00065             "where infile is an AFNI *.1D file (ASCII list of numbers arranged\n"
00066             "in columns); outfile will be a similar file, with each column being\n"
00067             "L2 normalized.\n"
00068            ) ;
00069       exit(0) ;
00070    }
00071 
00072    machdep() ;
00073 
00074    nopt = 1 ;
00075    nwin = strtol( argv[nopt] , NULL , 10 ) ;
00076    if( nwin < 16 ){
00077       fprintf(stderr,"** Illegal nwin!\n"); exit(1);
00078    }
00079 
00080    nopt++ ;
00081    if( nopt+1 >= argc ){
00082       fprintf(stderr,"** Need input and output filenames!\n");exit(1);
00083    }
00084 
00085    if( !THD_filename_ok(argv[nopt+1]) ){
00086       fprintf(stderr,"** Illegal output filename!\n"); exit(1);
00087    }
00088    if( THD_is_file(argv[nopt+1]) ){
00089       fprintf(stderr,"** Output file already exists!\n"); exit(1);
00090    }
00091 
00092    /* read input file */
00093 
00094    inim = mri_read_1D( argv[nopt] ) ;
00095    if( inim == NULL ){
00096       fprintf(stderr,"** Can't read input file!\n"); exit(1);
00097    }
00098 
00099    outim = ts_to_ftime( nwin , inim ) ;
00100    if( outim == NULL ){
00101       fprintf(stderr,"** Can't compute output!\n"); exit(1);
00102    }
00103 
00104    mri_write_1D( argv[nopt+1] , outim ) ;
00105    exit(0) ;
00106 }