Skip to content

AFNI/NIfTI Server

Sections
Personal tools
You are here: Home » AFNI » Documentation

Doxygen Source Code Documentation


Main Page   Alphabetical List   Data Structures   File List   Data Fields   Globals   Search  

ibeta.c File Reference

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

Go to the source code of this file.


Defines

#define NJ   9
#define NPT   ((1<<NJ)*10)
#define AITKEN(x, y, z)   ((x) - ((y)-(x))*((y)-(x)) / (((z)-(y))-((y)-(x))) )
#define DMP(tbl, jt)

Functions

void ibeta (double xcut, double a, double b, double *f1, double *flogx, double *flog1x)
int main (int argc, char *argv[])

Define Documentation

#define AITKEN x,
y,
z       ((x) - ((y)-(x))*((y)-(x)) / (((z)-(y))-((y)-(x))) )
 

#define DMP tbl,
jt   
 

Value:

do{ fprintf(stderr," table " #tbl ":") ;                              \
     for( jj=0 ; jj <= jt ; jj++ ) fprintf(stderr," %12.5g",tbl[jj]) ; \
     fprintf(stderr,"\n") ; } while(0)

#define NJ   9
 

Definition at line 6 of file ibeta.c.

Referenced by ibeta().

#define NPT   ((1<<NJ)*10)
 

Definition at line 9 of file ibeta.c.

Referenced by ibeta().


Function Documentation

void ibeta double    xcut,
double    a,
double    b,
double *    f1,
double *    flogx,
double *    flog1x
 

Definition at line 11 of file ibeta.c.

References a, NJ, and NPT.

Referenced by main(), and spmpar().

00013 {
00014    int nn , ii , jj , istep ;
00015    double dx , s1,sx,s1x , x, a1=a-1.0,b1=b-1.0 ;
00016    double as1[NJ+1] , asx[NJ+1] , as1x[NJ+1] ;
00017    double bs1[NJ-1] , bsx[NJ-1] , bs1x[NJ-1] ;
00018    double cs1[NJ-3] , csx[NJ-3] , cs1x[NJ-3] ;
00019    double xab[NPT] , xabx[NPT] , xab1x[NPT] ;
00020 
00021    if( xcut <= 0.0 || xcut >= 1.0 || a <= 0.0 || b <= 0.0 ) return ;
00022 
00023    if( f1 == NULL && flogx == NULL && flog1x == NULL ) return ;
00024 
00025    dx = xcut / NPT ;
00026    for( ii=1 ; ii <= NPT ; ii++ ){
00027       x = ii*dx ;
00028 /* fprintf(stderr,"ii=%d dx=%g a1=%g b1=%g x=%g\n",ii,dx,a1,b1,x) ; */
00029       xab[ii]   = pow(x,a1) * pow(1.0-x,b1) ;
00030       xabx[ii]  = xab[ii] * log(x) ;
00031       xab1x[ii] = xab[ii] * log(1.0-x) ;
00032 
00033 /* fprintf(stderr,"ii=%d x=%g xab=%g xabx=%g xab1x=%g a1=%g b1=%g dx=%g\n",
00034         ii,x,xab[ii],xabx[ii],xab1x[ii],a1,b1,dx) ; */
00035    }
00036 
00037    for( nn=NPT,istep=1,jj=NJ ; jj >= 0 ; jj--,istep*=2,nn/=2 ){
00038       s1 = sx = s1x = 0.0 ;
00039       for( ii=istep ; ii <= NPT ; ii+=istep ){
00040         s1 += xab[ii] ; sx += xabx[ii] ; s1x += xab1x[ii] ;
00041       }
00042       dx = xcut / nn ;
00043       as1[jj] = s1*dx ; asx[jj] = sx*dx ; as1x[jj] = s1x*dx ;
00044    }
00045 
00046 #undef  AITKEN
00047 #define AITKEN(x,y,z) ((x) - ((y)-(x))*((y)-(x)) / (((z)-(y))-((y)-(x))) )
00048 
00049    for( jj=0 ; jj <= NJ-2 ; jj++ ){
00050       bs1[jj]  = AITKEN( as1[jj]  , as1[jj+1]  , as1[jj+2]  ) ;
00051       bsx[jj]  = AITKEN( asx[jj]  , asx[jj+1]  , asx[jj+2]  ) ;
00052       bs1x[jj] = AITKEN( as1x[jj] , as1x[jj+1] , as1x[jj+2] ) ;
00053    }
00054 
00055    for( jj=0 ; jj <= NJ-4 ; jj++ ){
00056       cs1[jj]  = AITKEN( bs1[jj]  , bs1[jj+1]  , bs1[jj+2]  ) ;
00057       csx[jj]  = AITKEN( bsx[jj]  , bsx[jj+1]  , bsx[jj+2]  ) ;
00058       cs1x[jj] = AITKEN( bs1x[jj] , bs1x[jj+1] , bs1x[jj+2] ) ;
00059    }
00060 
00061 #define DMP(tbl,jt)                                                    \
00062  do{ fprintf(stderr," table " #tbl ":") ;                              \
00063      for( jj=0 ; jj <= jt ; jj++ ) fprintf(stderr," %12.5g",tbl[jj]) ; \
00064      fprintf(stderr,"\n") ; } while(0)
00065 
00066    fprintf(stderr,"ibeta(xc=%g,a=%g,b=%g)\n",xcut,a,b) ;
00067    DMP(as1,NJ) ; DMP(bs1,NJ-2) ; DMP(cs1,NJ-4) ;
00068    DMP(asx,NJ) ; DMP(bsx,NJ-2) ; DMP(csx,NJ-4) ;
00069    DMP(as1x,NJ) ; DMP(bs1x,NJ-2) ; DMP(cs1x,NJ-4) ;
00070 
00071    if( f1     != NULL ) *f1     = cs1[NJ-4]  ;
00072    if( flogx  != NULL ) *flogx  = csx[NJ-4]  ;
00073    if( flog1x != NULL ) *flog1x = cs1x[NJ-4] ;
00074 
00075    return ;
00076 }

int main int    argc,
char *    argv[]
 

\** File : SUMA.c

Author:
: Ziad Saad Date : Thu Dec 27 16:21:01 EST 2001
Purpose :

Input paramters :

Parameters:
param  Usage : SUMA ( )
Returns :
Returns:
Support :
See also:
OpenGL prog. Guide 3rd edition , varray.c from book's sample code
Side effects :

Definition at line 78 of file ibeta.c.

References a, argc, ibeta(), and strtod().

00079 {
00080    double xcut , a , b , f1,f2,f3 ;
00081 
00082    if( argc < 4 ) exit(0) ;
00083    xcut = strtod(argv[1],NULL) ;
00084    a    = strtod(argv[2],NULL) ;
00085    b    = strtod(argv[3],NULL) ;
00086 
00087    ibeta( xcut,a,b , &f1,&f2,&f3 ) ; exit(0) ;
00088 }
 

Powered by Plone

This site conforms to the following standards: