Doxygen Source Code Documentation
3dRegAna.c File Reference
#include <stdio.h>
#include <math.h>
#include "mrilib.h"
#include "matrix.h"
#include "RegAna.c"
Go to the source code of this file.
Data Structures | |
struct | model |
struct | RA_options |
Defines | |
#define | PROGRAM_NAME "3dRegAna" |
#define | PROGRAM_AUTHOR "B. Douglas Ward" |
#define | PROGRAM_INITIAL "10 Oct 1997" |
#define | PROGRAM_LATEST "02 Dec 2002" |
#define | SUFFIX ".3dregana" |
#define | DF "df -k ." |
#define | MAX_XVARS 101 |
#define | MAX_OBSERVATIONS 1000 |
#define | MAX_NAME_LENGTH THD_MAX_NAME |
#define | MEGA 1048576 |
#define | DOPEN(ds, name) |
#define | SUB_POINTER(ds, vv, ind, ptr) |
Typedefs | |
typedef model | model |
typedef RA_options | RA_options |
Functions | |
void | RA_error (char *message) |
void | display_help_menu () |
void | get_dimensions (RA_options *option_data) |
void | read_afni_data (RA_options *option_data, char *filename, int piece_len, int fim_offset, float *ffim) |
void | check_piece (char *filename) |
void | read_piece (char *filename, float *fin, int size) |
void | write_piece (char *filename, float *fout, int size) |
void | delete_piece (char *filename) |
void | allocate_pieces (matrix xdata, model *regmodel, RA_options *option_data, float ***yfimar, float **freg_piece, float **rsqr_piece, float ***coef_piece, float ***tcoef_piece) |
void | save_pieces (int piece, int piece_len, float *freg_piece, float *rsqr_piece, float **coef_piece, float **tcoef_piece) |
void | deallocate_pieces (int n, float ***yfimar, float **freg_piece, float **rsqr_piece, float ***coef_piece, float ***tcoef_piece) |
void | check_volume (char *filename, int num_pieces) |
void | read_volume (char *filename, float *volume, int nxyz, int piece_size, int num_pieces) |
void | delete_volume (char *filename, int nxyz, int piece_size, int num_pieces) |
void | initialize_options (model *regmodel, RA_options *option_data) |
void | sort_model_indices (model *regmodel) |
void | get_inputs (int argc, char **argv, matrix *xdata, model *regmodel, RA_options *option_data) |
void | count_volumes_and_files (model *regmodel, RA_options *option_data) |
void | break_into_pieces (int num_datasets, RA_options *option_data) |
void | identify_repeats (matrix *xdata, model *regmodel, RA_options *option_data) |
void | check_for_valid_inputs (matrix *xdata, model *regmodel, RA_options *option_data) |
void | check_one_output_file (RA_options *option_data, char *filename) |
void | check_output_files (RA_options *option_data) |
void | check_temporary_files (matrix *xdata, model *regmodel, RA_options *option_data) |
void | check_disk_space (RA_options *option_data) |
void | initialize_program (int argc, char **argv, matrix *xdata, model *regmodel, RA_options *option_data) |
void | init_regression_analysis (int p, int q, int *flist, int *rlist, matrix xdata, matrix *x_full, matrix *xtxinv_full, matrix *xtxinvxt_full, matrix *x_base, matrix *xtxinvxt_base, matrix *x_rdcd, matrix *xtxinvxt_rdcd) |
void | regression_analysis (int N, int p, int q, matrix x_full, matrix xtxinv_full, matrix xtxinvxt_full, matrix x_base, matrix xtxinvxt_base, matrix x_rdcd, matrix xtxinvxt_rdcd, vector y, float rms_min, int *levels, int *counts, int c, float flofmax, float *flof, vector *coef_full, vector *scoef_full, vector *tcoef_full, float *freg, float *rsqr) |
void | save_voxel (int iv, vector y, float fdisp, model *regmodel, float flof, vector coef, vector tcoef, float freg, float rsqr, float *freg_piece, float *rsqr_piece, float **coef_piece, float **tcoef_piece) |
void | calculate_results (matrix xdata, model *regmodel, RA_options *option_data) |
float | EDIT_coerce_autoscale_new (int nxyz, int itype, void *ivol, int otype, void *ovol) |
void | write_afni_data (RA_options *option_data, char *filename, float *ffim, float *ftr, int func_type, int numdof, int dendof) |
void | write_bucket_data (matrix xdata, model *regmodel, RA_options *option_data) |
void | output_results (matrix xdata, model *regmodel, RA_options *option_data) |
void | terminate_program (matrix *xdata, model *regmodel, RA_options *option_data) |
int | main (int argc, char **argv) |
Variables | |
char * | commandline = NULL |
Define Documentation
|
Definition at line 104 of file 3dRegAna.c. Referenced by check_disk_space(). |
|
Value: do{ int pv ; (ds) = THD_open_dataset((name)) ; \ if( !ISVALID_3DIM_DATASET((ds)) ){ \ fprintf(stderr,"*** Can't open dataset: %s\n",(name)) ; exit(1) ; } \ if( (ds)->daxes->nxx!=nx || (ds)->daxes->nyy!=ny || \ (ds)->daxes->nzz!=nz ){ \ fprintf(stderr,"*** Axes mismatch: %s\n",(name)) ; exit(1) ; } \ if( DSET_NVALS((ds)) != 1 ){ \ fprintf(stderr,"*** Must specify 1 sub-brick: %s\n",(name));exit(1);}\ THD_load_datablock( (ds)->dblk ) ; \ pv = DSET_PRINCIPAL_VALUE((ds)) ; \ if( DSET_ARRAY((ds),pv) == NULL ){ \ fprintf(stderr,"*** Can't access data: %s\n",(name)) ; exit(1); } \ if( DSET_BRICK_TYPE((ds),pv) == MRI_complex ){ \ fprintf(stderr,"*** Can't use complex data: %s\n",(name)); exit(1); \ } \ break ; } while (0) Definition at line 269 of file 3dRegAna.c. |
|
Definition at line 117 of file 3dRegAna.c. Referenced by check_disk_space(), check_piece(), check_temporary_files(), check_volume(), delete_piece(), delete_volume(), get_inputs(), output_results(), read_piece(), read_volume(), save_pieces(), terminate_program(), write_bucket_data(), and write_piece(). |
|
Definition at line 116 of file 3dRegAna.c. Referenced by get_inputs(), initialize_options(), and terminate_program(). |
|
Definition at line 115 of file 3dRegAna.c. Referenced by allocate_pieces(), check_output_files(), deallocate_pieces(), get_inputs(), init_regression_analysis(), initialize_options(), save_pieces(), and terminate_program(). |
|
Definition at line 118 of file 3dRegAna.c. Referenced by break_into_pieces(). |
|
Definition at line 60 of file 3dRegAna.c. Referenced by main(). |
|
Definition at line 61 of file 3dRegAna.c. Referenced by main(). |
|
Definition at line 62 of file 3dRegAna.c. Referenced by main(). |
|
Definition at line 59 of file 3dRegAna.c. Referenced by get_inputs(), initialize_program(), main(), and RA_error(). |
|
Value: do{ switch( DSET_BRICK_TYPE((ds),(vv)) ){ \ default: fprintf(stderr,"\n*** Illegal datum! ***\n");exit(1); \ case MRI_short:{ short * fim = (short *) DSET_ARRAY((ds),(vv)) ; \ (ptr) = (void *)( fim + (ind) ) ; \ } break ; \ case MRI_byte:{ byte * fim = (byte *) DSET_ARRAY((ds),(vv)) ; \ (ptr) = (void *)( fim + (ind) ) ; \ } break ; \ case MRI_float:{ float * fim = (float *) DSET_ARRAY((ds),(vv)) ; \ (ptr) = (void *)( fim + (ind) ) ; \ } break ; } break ; } while(0) Definition at line 292 of file 3dRegAna.c. |
|
Definition at line 66 of file 3dRegAna.c. Referenced by check_piece(), delete_piece(), read_piece(), and write_piece(). |
Typedef Documentation
|
|
|
|
Function Documentation
|
Definition at line 537 of file 3dRegAna.c. References FUNC_FIM_TYPE, FUNC_THR_TYPE, i, malloc, MAX_XVARS, MTEST, and p. Referenced by calculate_results().
00549 { 00550 int n; /* number of datasets */ 00551 int p; /* number of parameters in the full model */ 00552 int * flist = NULL; /* list of parameters in the full model */ 00553 int i; /* dataset index */ 00554 int ip; /* parameter index */ 00555 int ix; /* x-variable index */ 00556 int piece_size; /* number of voxels in datasset piece */ 00557 int ib; /* sub-brick index */ 00558 int nbricks; /* number of sub-bricks in bucket dataset */ 00559 00560 00561 /*----- initialize local variables -----*/ 00562 n = xdata.rows; 00563 p = regmodel->p; 00564 flist = regmodel->flist; 00565 piece_size = option_data->piece_size; 00566 nbricks = option_data->numbricks; 00567 00568 00569 /*----- allocate memory space for input data -----*/ 00570 *yfimar = (float **) malloc (sizeof(float *) * n); 00571 MTEST (*yfimar); 00572 for (i = 0; i < n; i++) 00573 { 00574 (*yfimar)[i] = (float *) malloc(sizeof(float) * piece_size); 00575 MTEST ((*yfimar)[i]); 00576 } 00577 00578 00579 /*----- allocate memory space for F-statistics data -----*/ 00580 for (ip = 0; ip < p; ip++) 00581 { 00582 ix = flist[ip]; 00583 00584 if (option_data->fcoef_filename[ix] != NULL) 00585 { 00586 if (*freg_piece == NULL) 00587 { 00588 *freg_piece = (float *) malloc (sizeof(float) * piece_size); 00589 MTEST (*freg_piece); 00590 } 00591 } 00592 } 00593 00594 for (ib = 0; ib < nbricks; ib++) 00595 { 00596 if (option_data->brick_type[ib] == FUNC_FT_TYPE) 00597 { 00598 if (*freg_piece == NULL) 00599 { 00600 *freg_piece = (float *) malloc (sizeof(float) * piece_size); 00601 MTEST (*freg_piece); 00602 } 00603 } 00604 } 00605 00606 00607 /*----- allocate memory space for coef. of mult. determination R^2 -----*/ 00608 for (ip = 0; ip < p; ip++) 00609 { 00610 ix = flist[ip]; 00611 00612 if (option_data->rcoef_filename[ix] != NULL) 00613 { 00614 if (*rsqr_piece == NULL) 00615 { 00616 *rsqr_piece = (float *) malloc (sizeof(float) * piece_size); 00617 MTEST (*rsqr_piece); 00618 } 00619 } 00620 } 00621 00622 for (ib = 0; ib < nbricks; ib++) 00623 { 00624 if (option_data->brick_type[ib] == FUNC_THR_TYPE) 00625 { 00626 if (*rsqr_piece == NULL) 00627 { 00628 *rsqr_piece = (float *) malloc (sizeof(float) * piece_size); 00629 MTEST (*rsqr_piece); 00630 } 00631 } 00632 } 00633 00634 00635 /*----- allocate memory space for parameters and t-stats -----*/ 00636 *coef_piece = (float **) malloc (sizeof(float *) * MAX_XVARS); 00637 MTEST (*coef_piece); 00638 00639 *tcoef_piece = (float **) malloc (sizeof(float *) * MAX_XVARS); 00640 MTEST (*tcoef_piece); 00641 00642 for (ix = 0; ix < MAX_XVARS; ix++) 00643 { 00644 (*coef_piece)[ix] = NULL; 00645 (*tcoef_piece)[ix] = NULL; 00646 } 00647 00648 for (ip = 0; ip < p; ip++) 00649 { 00650 ix = flist[ip]; 00651 00652 if ( (option_data->fcoef_filename[ix] != NULL) 00653 || (option_data->rcoef_filename[ix] != NULL) 00654 || (option_data->tcoef_filename[ix] != NULL)) 00655 { 00656 (*coef_piece)[ix] = 00657 (float *) malloc (sizeof(float) * piece_size); 00658 MTEST ((*coef_piece)[ix]); 00659 } 00660 00661 if (option_data->tcoef_filename[ix] != NULL) 00662 { 00663 (*tcoef_piece)[ix] = 00664 (float *) malloc (sizeof(float) * piece_size); 00665 MTEST ((*tcoef_piece)[ix]); 00666 } 00667 } 00668 00669 for (ib = 0; ib < nbricks; ib++) 00670 { 00671 if (option_data->brick_type[ib] == FUNC_FIM_TYPE) 00672 { 00673 ix = option_data->brick_coef[ib]; 00674 if ((*coef_piece)[ix] == NULL) 00675 { 00676 (*coef_piece)[ix] = (float *) malloc (sizeof(float)*piece_size); 00677 MTEST ((*coef_piece)[ix]); 00678 } 00679 } 00680 } 00681 00682 for (ib = 0; ib < nbricks; ib++) 00683 { 00684 if (option_data->brick_type[ib] == FUNC_TT_TYPE) 00685 { 00686 ix = option_data->brick_coef[ib]; 00687 if ((*tcoef_piece)[ix] == NULL) 00688 { 00689 (*tcoef_piece)[ix] = (float *) malloc (sizeof(float)*piece_size); 00690 MTEST ((*tcoef_piece)[ix]); 00691 } 00692 } 00693 } 00694 00695 } |
|
Definition at line 1636 of file 3dRegAna.c. References MEGA. Referenced by initialize_program().
01641 { 01642 int num_vols; /* number of output volumes */ 01643 01644 01645 /*----- count number of distinct floating point volumes required -----*/ 01646 num_vols = option_data->numf + option_data->numr 01647 + option_data->numt + option_data->numc; 01648 01649 /*----- calculate number of voxels per piece -----*/ 01650 option_data->piece_size = option_data->workmem * MEGA 01651 / ((num_datasets + num_vols) * sizeof(float)); 01652 if (option_data->piece_size > option_data->nxyz) 01653 option_data->piece_size = option_data->nxyz; 01654 01655 /*----- calculate number of pieces per dataset -----*/ 01656 option_data->num_pieces = (option_data->nxyz + option_data->piece_size - 1) 01657 / option_data->piece_size; 01658 01659 printf ("num_pieces = %d piece_size = %d \n", 01660 option_data->num_pieces, option_data->piece_size); 01661 01662 } |
|
Definition at line 2321 of file 3dRegAna.c. References allocate_pieces(), deallocate_pieces(), vector::elts, i, init_regression_analysis(), matrix_destroy(), matrix_initialize(), matrix_print(), p, q, read_afni_data(), regression_analysis(), save_pieces(), save_voxel(), vector_create(), vector_destroy(), and vector_initialize().
02327 { 02328 int * flist = NULL; /* list of parameters in the full model */ 02329 int * rlist = NULL; /* list of parameters in the rdcd model */ 02330 int n; /* number of data points */ 02331 int p; /* number of parameters in the full model */ 02332 int q; /* number of parameters in the rdcd model */ 02333 float flof; /* F-statistic for lack of fit */ 02334 float freg; /* F-statistic for the full regression model */ 02335 float rsqr; /* coeff. of multiple determination R^2 */ 02336 vector coef; /* regression parameters */ 02337 vector scoef; /* std. devs. for regression parameters */ 02338 vector tcoef; /* t-statistics for regression parameters */ 02339 02340 matrix x_full; /* extracted X matrix for full model */ 02341 matrix xtxinv_full; /* matrix: 1/(X'X) for full model */ 02342 matrix xtxinvxt_full; /* matrix: (1/(X'X))X' for full model */ 02343 matrix x_base; /* extracted X matrix for baseline model */ 02344 matrix xtxinvxt_base; /* matrix: (1/(X'X))X' for baseline model */ 02345 matrix x_rdcd; /* extracted X matrix for reduced model */ 02346 matrix xtxinvxt_rdcd; /* matrix: (1/(X'X))X' for reduced model */ 02347 vector y; /* vector of measured data */ 02348 02349 int i; /* dataset index */ 02350 int nxyz; /* number of voxels per dataset */ 02351 int num_datasets; /* total number of datasets */ 02352 int piece_size; /* number of voxels in dataset piece */ 02353 int num_pieces; /* dataset is divided into this many pieces */ 02354 int piece; /* piece index */ 02355 int piece_len; /* number of voxels in current piece */ 02356 int fim_offset; /* array offset to current piece */ 02357 int ivox; /* index of voxel in current piece */ 02358 int nvox; /* index of voxel within entire volume */ 02359 02360 float ** yfimar = NULL; /* array of pieces of Y-datasets */ 02361 float * freg_piece = NULL; /* piece for F-statistics */ 02362 float * rsqr_piece = NULL; /* piece for R^2 */ 02363 float ** coef_piece = NULL; /* pieces for regression coefficients */ 02364 float ** tcoef_piece = NULL; /* pieces for t-statistics */ 02365 02366 02367 /*----- Initialize matrices and vectors -----*/ 02368 matrix_initialize (&x_full); 02369 matrix_initialize (&xtxinv_full); 02370 matrix_initialize (&xtxinvxt_full); 02371 matrix_initialize (&x_base); 02372 matrix_initialize (&xtxinvxt_base); 02373 matrix_initialize (&x_rdcd); 02374 matrix_initialize (&xtxinvxt_rdcd); 02375 vector_initialize (&coef); 02376 vector_initialize (&scoef); 02377 vector_initialize (&tcoef); 02378 vector_initialize (&y); 02379 02380 02381 /*----- initialize local variables -----*/ 02382 n = xdata.rows; 02383 p = regmodel->p; 02384 flist = regmodel->flist; 02385 q = regmodel->q; 02386 rlist = regmodel->rlist; 02387 nxyz = option_data->nxyz; 02388 piece_size = option_data->piece_size; 02389 num_pieces = option_data->num_pieces; 02390 02391 02392 /*----- allocate memory space for pieces -----*/ 02393 allocate_pieces (xdata, regmodel, option_data, &yfimar, 02394 &freg_piece, &rsqr_piece, &coef_piece, &tcoef_piece); 02395 02396 02397 /*----- initialization for the regression analysis -----*/ 02398 init_regression_analysis (p, q, flist, rlist, xdata, 02399 &x_full, &xtxinv_full, &xtxinvxt_full, 02400 &x_base, &xtxinvxt_base, &x_rdcd, &xtxinvxt_rdcd); 02401 02402 02403 vector_create (n, &y); 02404 02405 02406 if (option_data->fdisp >= 0) 02407 { 02408 printf ("\n"); 02409 printf ("X matrix: \n"); 02410 matrix_print (xdata); 02411 } 02412 02413 02414 /*----- loop over the pieces of the input datasets -----*/ 02415 nvox = -1; 02416 for (piece = 0; piece < num_pieces; piece++) 02417 { 02418 printf ("piece = %d \n", piece); 02419 fim_offset = piece * piece_size; 02420 if (piece < num_pieces-1) 02421 piece_len = piece_size; 02422 else 02423 piece_len = nxyz - fim_offset; 02424 02425 /*----- read in the Y-data -----*/ 02426 for (i = 0; i < n; i++) 02427 read_afni_data (option_data, option_data->yname[i], 02428 piece_len, fim_offset, yfimar[i]); 02429 02430 02431 /*----- loop over voxels in this piece -----*/ 02432 for (ivox = 0; ivox < piece_len; ivox++) 02433 { 02434 nvox += 1; 02435 02436 02437 /*----- extract Y-data for this voxel -----*/ 02438 for (i = 0; i < n; i++) 02439 y.elts[i] = yfimar[i][ivox]; 02440 02441 02442 /*----- calculate results for this voxel -----*/ 02443 regression_analysis (n, p, q, 02444 x_full, xtxinv_full, xtxinvxt_full, x_base, 02445 xtxinvxt_base, x_rdcd, xtxinvxt_rdcd, 02446 y, option_data->rms_min, option_data->levels, 02447 option_data->counts, option_data->c, 02448 option_data->flofmax, &flof, 02449 &coef, &scoef, &tcoef, &freg, &rsqr); 02450 02451 02452 /*----- save results for this voxel -----*/ 02453 save_voxel (ivox, y, option_data->fdisp, 02454 regmodel, flof, coef, tcoef, freg, rsqr, 02455 freg_piece, rsqr_piece, coef_piece, tcoef_piece); 02456 02457 02458 } /* loop over voxels within this piece */ 02459 02460 02461 /*----- save piece data into external files -----*/ 02462 save_pieces (piece, piece_len, 02463 freg_piece, rsqr_piece, coef_piece, tcoef_piece); 02464 02465 02466 } /* loop over pieces */ 02467 02468 02469 /*----- deallocate memory for pieces -----*/ 02470 deallocate_pieces (n, &yfimar, &freg_piece, &rsqr_piece, 02471 &coef_piece, &tcoef_piece); 02472 02473 02474 /*----- Dispose of matrices -----*/ 02475 vector_destroy (&y); 02476 vector_destroy (&tcoef); 02477 vector_destroy (&scoef); 02478 vector_destroy (&coef); 02479 matrix_destroy (&xtxinvxt_rdcd); 02480 matrix_destroy (&x_rdcd); 02481 matrix_destroy (&xtxinvxt_base); 02482 matrix_destroy (&x_base); 02483 matrix_destroy (&xtxinvxt_full); 02484 matrix_destroy (&xtxinv_full); 02485 matrix_destroy (&x_full); 02486 02487 } |
|
Definition at line 1973 of file 3dRegAna.c. References DF, and MAX_NAME_LENGTH. Referenced by initialize(), and initialize_program().
01977 { 01978 char ch; /* user response */ 01979 int nxyz; /* number of voxels per image */ 01980 int nmax; /* maximum number of disk files */ 01981 char filename[MAX_NAME_LENGTH]; /* output file name */ 01982 01983 01984 /*----- initialize local variables -----*/ 01985 nxyz = option_data->nxyz; 01986 01987 /*----- first, determine space required for temporary volume data -----*/ 01988 nmax = 4 * nxyz * (option_data->numf + option_data->numr + option_data->numt 01989 + option_data->numc); 01990 01991 /*----- determine additional space required for output files -----*/ 01992 nmax += 4 * nxyz * option_data->numfiles; 01993 01994 printf("\nThis problem requires approximately %d MB of free disk space.\n", 01995 (nmax/1000000) + 1); 01996 01997 01998 /*----- Determine how much disk space is available. -----*/ 01999 printf ("Summary of available disk space: \n\n"); 02000 system (DF); 02001 printf ("\nDo you wish to proceed? "); 02002 ch = getchar(); 02003 if ((ch != 'Y') && (ch != 'y')) exit(0); 02004 printf ("\n"); 02005 } |
|
Definition at line 1755 of file 3dRegAna.c. References RA_error.
01761 { 01762 int ib; /* sub-brick index */ 01763 01764 01765 /*----- check data set dimensions -----*/ 01766 if (xdata->cols < 2) 01767 RA_error ("Too few X variables "); 01768 if (xdata->rows < 3) 01769 RA_error ("Too few data sets for Y-observations "); 01770 01771 /*----- check model dimensions -----*/ 01772 if (regmodel->q < 1) 01773 RA_error ("Reduced regression model is too small"); 01774 if (regmodel->p <= regmodel->q) 01775 RA_error ("Full regression model is too small"); 01776 if (xdata->rows <= regmodel->p) 01777 RA_error ("Too few data sets for fitting full regression model "); 01778 if (regmodel->rlist[0] != 0) 01779 RA_error ("Reduced model must include variable index 0"); 01780 01781 01782 /*----- check for repeat observations -----*/ 01783 if (option_data->flofmax >= 0.0) 01784 { 01785 if (xdata->rows <= option_data->c) 01786 RA_error ("Cannot conduct F-test for lack of fit (no repeat obs.) "); 01787 if (option_data->c <= regmodel->p) 01788 RA_error ("Cannot conduct F-test for lack of fit (c <= p) "); 01789 } 01790 01791 /*----- check bucket dataset parameters -----*/ 01792 if (option_data->numbricks > 0) 01793 for (ib = 0; ib < option_data->numbricks; ib++) 01794 { 01795 if (option_data->brick_type[ib] < 0) 01796 RA_error ("Must specify contents of each brick in the bucket"); 01797 } 01798 } |
|
Definition at line 1807 of file 3dRegAna.c. References ADN_directory_name, ADN_label1, ADN_none, ADN_prefix, ADN_self_name, ADN_type, THD_3dim_dataset::dblk, THD_datablock::diskptr, EDIT_dset_items(), EDIT_empty_copy(), GEN_FUNC_TYPE, HEAD_FUNC_TYPE, THD_diskptr::header_name, ISHEAD, ISVALID_3DIM_DATASET, THD_delete_3dim_dataset(), THD_is_file(), and THD_open_dataset().
01812 { 01813 THD_3dim_dataset * dset=NULL; /* input afni data set pointer */ 01814 THD_3dim_dataset * new_dset=NULL; /* output afni data set pointer */ 01815 int ierror; /* number of errors in editing data */ 01816 01817 01818 /*----- read first dataset -----*/ 01819 dset = THD_open_dataset (option_data->first_dataset ) ; 01820 if( ! ISVALID_3DIM_DATASET(dset) ){ 01821 fprintf(stderr,"*** Unable to open dataset file %s\n", 01822 option_data->first_dataset); 01823 exit(1) ; 01824 } 01825 01826 /*-- make an empty copy of this dataset, for eventual output --*/ 01827 new_dset = EDIT_empty_copy( dset ) ; 01828 01829 01830 ierror = EDIT_dset_items( new_dset , 01831 ADN_prefix , filename , 01832 ADN_label1 , filename , 01833 ADN_directory_name , option_data->session , 01834 ADN_self_name , filename , 01835 ADN_type , ISHEAD(dset) ? HEAD_FUNC_TYPE : 01836 GEN_FUNC_TYPE , 01837 ADN_none ) ; 01838 01839 if( ierror > 0 ){ 01840 fprintf(stderr, 01841 "*** %d errors in attempting to create output dataset!\n", ierror ) ; 01842 exit(1) ; 01843 } 01844 01845 if( THD_is_file(new_dset->dblk->diskptr->header_name) ){ 01846 fprintf(stderr, 01847 "*** Output dataset file %s already exists--cannot continue!\a\n", 01848 new_dset->dblk->diskptr->header_name ) ; 01849 exit(1) ; 01850 } 01851 01852 /*----- deallocate memory -----*/ 01853 THD_delete_3dim_dataset( dset , False ) ; dset = NULL ; 01854 THD_delete_3dim_dataset( new_dset , False ) ; new_dset = NULL ; 01855 01856 } |
|
Definition at line 1865 of file 3dRegAna.c. References check_one_output_file(), and MAX_XVARS.
01869 { 01870 int ix; /* x-variable index */ 01871 01872 01873 for (ix = 0; ix < MAX_XVARS; ix++) 01874 { 01875 if (option_data->fcoef_filename[ix] != NULL) 01876 check_one_output_file (option_data, option_data->fcoef_filename[ix]); 01877 if (option_data->rcoef_filename[ix] != NULL) 01878 check_one_output_file (option_data, option_data->rcoef_filename[ix]); 01879 if (option_data->tcoef_filename[ix] != NULL) 01880 check_one_output_file (option_data, option_data->tcoef_filename[ix]); 01881 } 01882 01883 if (option_data->bucket_filename != NULL) 01884 check_one_output_file (option_data, option_data->bucket_filename); 01885 01886 } |
|
Definition at line 387 of file 3dRegAna.c. References far, MAX_NAME_LENGTH, RA_error, and SUFFIX. Referenced by check_volume().
00391 { 00392 FILE * far = NULL; /* pointer to temporary file */ 00393 char sfilename[MAX_NAME_LENGTH]; /* name of temporary file */ 00394 char message[MAX_NAME_LENGTH]; /* error message */ 00395 00396 00397 /*----- see if piece file already exists -----*/ 00398 strcpy (sfilename, filename); 00399 strcat (sfilename, SUFFIX); 00400 far = fopen (sfilename, "r"); 00401 if (far != NULL) 00402 { 00403 sprintf (message, "temporary file %s already exists. ", sfilename); 00404 RA_error (message); 00405 } 00406 00407 } |
|
Definition at line 1895 of file 3dRegAna.c. References check_volume(), MAX_NAME_LENGTH, and p. Referenced by initialize(), and initialize_program().
01901 { 01902 int p; /* number of parameters in full model */ 01903 int ip, ix; /* parameter index */ 01904 int num_pieces; /* dataset is divided into this many pieces */ 01905 char filename[MAX_NAME_LENGTH]; /* name for temporary data file */ 01906 01907 01908 /*----- initialize local variables -----*/ 01909 p = regmodel->p; 01910 num_pieces = option_data->num_pieces; 01911 01912 01913 /*----- check for F-statistics data file -----*/ 01914 if (option_data->numf > 0) 01915 { 01916 strcpy (filename, "freg"); 01917 check_volume (filename, num_pieces); 01918 } 01919 01920 01921 /*----- check for R^2 data file -----*/ 01922 if (option_data->numr > 0) 01923 { 01924 strcpy (filename, "rsqr"); 01925 check_volume (filename, num_pieces); 01926 } 01927 01928 01929 /*----- check for t-statistics data files -----*/ 01930 if (option_data->numt > 0) 01931 { 01932 for (ip = 0; ip < p; ip++) 01933 { 01934 ix = regmodel->flist[ip]; 01935 01936 if (option_data->tcoef_filename[ix] != NULL) 01937 { 01938 sprintf (filename, "tcoef.%d", ix); 01939 check_volume (filename, num_pieces); 01940 } 01941 } 01942 } 01943 01944 01945 /*----- check for regression coefficients data files -----*/ 01946 if (option_data->numc > 0) 01947 { 01948 for (ip = 0; ip < p; ip++) 01949 { 01950 ix = regmodel->flist[ip]; 01951 01952 if ( (option_data->fcoef_filename[ix] != NULL) 01953 || (option_data->rcoef_filename[ix] != NULL) 01954 || (option_data->tcoef_filename[ix] != NULL) ) 01955 { 01956 sprintf (filename, "coef.%d", ix); 01957 check_volume (filename, num_pieces); 01958 } 01959 } 01960 } 01961 01962 01963 } |
|
Definition at line 841 of file 3dRegAna.c. References check_piece(), and MAX_NAME_LENGTH. Referenced by check_temporary_files().
00846 { 00847 int piece; /* piece index */ 00848 char sfilename[MAX_NAME_LENGTH]; /* name for temporary data file */ 00849 00850 00851 /*----- loop over the temporary data file pieces -----*/ 00852 for (piece = 0; piece < num_pieces; piece++) 00853 { 00854 /*----- check for piece data file -----*/ 00855 sprintf (sfilename, "%s.p%d", filename, piece); 00856 check_piece (sfilename); 00857 00858 } /* loop over pieces */ 00859 00860 } |
|
Definition at line 1539 of file 3dRegAna.c. References FUNC_FIM_TYPE, FUNC_THR_TYPE, and p. Referenced by initialize_program().
01544 { 01545 int ip; /* parameter index */ 01546 int ix; /* x-variable index */ 01547 int ib; /* sub-brick index */ 01548 int p; /* number of parameters in the model */ 01549 int nbricks; /* number of bricks in bucket dataset */ 01550 01551 01552 /*----- initialize local variables -----*/ 01553 nbricks = option_data->numbricks; 01554 p = regmodel->p; 01555 01556 01557 /*----- count number of volumes and files for F-statistics -----*/ 01558 for (ip = 0; ip < p; ip++) 01559 { 01560 ix = regmodel->flist[ip]; 01561 01562 if (option_data->fcoef_filename[ix] != NULL) 01563 { 01564 option_data->numf = 1; 01565 option_data->numfiles += 1; 01566 } 01567 } 01568 01569 for (ib = 0; ib < nbricks; ib++) 01570 if (option_data->brick_type[ib] == FUNC_FT_TYPE) 01571 option_data->numf = 1; 01572 01573 01574 /*----- count number of volumes and files for R^2 -----*/ 01575 for (ip = 0; ip < p; ip++) 01576 { 01577 ix = regmodel->flist[ip]; 01578 01579 if (option_data->rcoef_filename[ix] != NULL) 01580 { 01581 option_data->numr = 1; 01582 option_data->numfiles += 1; 01583 } 01584 } 01585 01586 for (ib = 0; ib < nbricks; ib++) 01587 if (option_data->brick_type[ib] == FUNC_THR_TYPE) 01588 option_data->numr = 1; 01589 01590 01591 /*----- count number of volumes and files for t-statistics -----*/ 01592 for (ip = 0; ip < p; ip++) 01593 { 01594 ix = regmodel->flist[ip]; 01595 01596 if (option_data->tcoef_filename[ix] != NULL) 01597 { 01598 option_data->numt += 1; 01599 option_data->numfiles += 1; 01600 } 01601 01602 else 01603 for (ib = 0; ib < nbricks; ib++) 01604 if ((option_data->brick_type[ib] == FUNC_TT_TYPE) 01605 && (option_data->brick_coef[ib] == ix)) 01606 option_data->numt += 1; 01607 } 01608 01609 01610 /*----- count number of volumes for regression coefficients -----*/ 01611 for (ip = 0; ip < p; ip++) 01612 { 01613 ix = regmodel->flist[ip]; 01614 01615 if ( (option_data->fcoef_filename[ix] != NULL) 01616 || (option_data->rcoef_filename[ix] != NULL) 01617 || (option_data->tcoef_filename[ix] != NULL)) 01618 option_data->numc += 1; 01619 01620 else 01621 for (ib = 0; ib < nbricks; ib++) 01622 if ((option_data->brick_type[ib] == FUNC_FIM_TYPE) 01623 && (option_data->brick_coef[ib] == ix)) 01624 option_data->numc += 1; 01625 } 01626 01627 } |
|
Definition at line 767 of file 3dRegAna.c. References free, i, and MAX_XVARS. Referenced by calculate_results().
00776 { 00777 int i; /* dataset index */ 00778 int ip; /* parameter index */ 00779 00780 00781 /*----- deallocate memory space for input data -----*/ 00782 if (*yfimar != NULL) 00783 { 00784 for (i = 0; i < n; i++) 00785 { 00786 free ((*yfimar)[i]); (*yfimar)[i] = NULL; 00787 } 00788 free (*yfimar); 00789 *yfimar = NULL; 00790 } 00791 00792 /*----- deallocate memory space for F-statistics data -----*/ 00793 if (*freg_piece != NULL) 00794 { 00795 free (*freg_piece); 00796 *freg_piece = NULL; 00797 } 00798 00799 /*----- deallocate memory space for coef. of mult. determination R^2 -----*/ 00800 if (*rsqr_piece != NULL) 00801 { 00802 free (*rsqr_piece); 00803 *rsqr_piece = NULL; 00804 } 00805 00806 /*----- deallocate memory space for regression coefficients -----*/ 00807 if (*coef_piece != NULL) 00808 { 00809 for (ip = 0; ip < MAX_XVARS; ip++) 00810 if ((*coef_piece)[ip] != NULL) 00811 { 00812 free ((*coef_piece)[ip]); 00813 (*coef_piece)[ip] = NULL; 00814 } 00815 free (*coef_piece); 00816 *coef_piece = NULL; 00817 } 00818 00819 /*----- deallocate memory space for t-statistics -----*/ 00820 if (*tcoef_piece != NULL) 00821 { 00822 for (ip = 0; ip < MAX_XVARS; ip++) 00823 if ((*tcoef_piece)[ip] != NULL) 00824 { 00825 free ((*tcoef_piece)[ip]); 00826 (*tcoef_piece)[ip] = NULL; 00827 } 00828 free (*tcoef_piece); 00829 *tcoef_piece = NULL; 00830 } 00831 00832 } |
|
Definition at line 514 of file 3dRegAna.c. References MAX_NAME_LENGTH, and SUFFIX. Referenced by delete_volume().
00518 { 00519 char sfilename[MAX_NAME_LENGTH]; /* file name */ 00520 00521 /*----- build file name -----*/ 00522 strcpy (sfilename, filename); 00523 strcat (sfilename, SUFFIX); 00524 00525 /*----- delete 3d data file -----*/ 00526 remove (sfilename); 00527 00528 } |
|
Definition at line 912 of file 3dRegAna.c. References delete_piece(), and MAX_NAME_LENGTH. Referenced by terminate_program(), and write_bucket_data().
00919 { 00920 int piece; /* piece index */ 00921 char sfilename[MAX_NAME_LENGTH]; /* name for temporary data file */ 00922 00923 00924 /*----- loop over the temporary data file pieces -----*/ 00925 for (piece = 0; piece < num_pieces; piece++) 00926 { 00927 00928 /*----- delete the piece data file -----*/ 00929 sprintf (sfilename, "%s.p%d", filename, piece); 00930 delete_piece (sfilename); 00931 00932 } /* loop over pieces */ 00933 00934 } |
|
Definition at line 194 of file 3dRegAna.c. References MASTER_SHORTHELP_STRING.
00195 { 00196 printf 00197 ( 00198 "This program performs multiple linear regression analysis. \n\n" 00199 "Usage: \n" 00200 "3dRegAna \n" 00201 "-rows n number of input datasets \n" 00202 "-cols m number of X variables \n" 00203 "-xydata X11 X12 ... X1m filename X variables and Y observations \n" 00204 " . . \n" 00205 " . . \n" 00206 " . . \n" 00207 "-xydata Xn1 Xn2 ... Xnm filename X variables and Y observations \n" 00208 " \n" 00209 "-model i1 ... iq : j1 ... jr definition of linear regression model; \n" 00210 " reduced model: \n" 00211 " Y = f(Xj1,...,Xjr) \n" 00212 " full model: \n" 00213 " Y = f(Xj1,...,Xjr,Xi1,...,Xiq) \n" 00214 " \n" 00215 "[-diskspace] print out disk space required for program execution\n" 00216 "[-workmem mega] number of megabytes of RAM to use for statistical \n" 00217 " workspace (default = 12) \n" 00218 "[-rmsmin r] r = minimum rms error to reject constant model \n" 00219 "[-fdisp fval] display (to screen) results for those voxels \n" 00220 " whose F-statistic is > fval \n" 00221 " \n" 00222 "[-flof alpha] alpha = minimum p value for F due to lack of fit \n" 00223 " \n" 00224 " \n" 00225 "The following commands generate individual AFNI 2 sub-brick datasets: \n" 00226 " \n" 00227 "[-fcoef k prefixname] estimate of kth regression coefficient \n" 00228 " along with F-test for the regression \n" 00229 " is written to AFNI `fift' dataset \n" 00230 "[-rcoef k prefixname] estimate of kth regression coefficient \n" 00231 " along with coef. of mult. deter. R^2 \n" 00232 " is written to AFNI `fith' dataset \n" 00233 "[-tcoef k prefixname] estimate of kth regression coefficient \n" 00234 " along with t-test for the coefficient \n" 00235 " is written to AFNI `fitt' dataset \n" 00236 " \n" 00237 " \n" 00238 "The following commands generate one AFNI 'bucket' type dataset: \n" 00239 " \n" 00240 "[-bucket n prefixname] create one AFNI 'bucket' dataset having \n" 00241 " n sub-bricks; n=0 creates default output;\n" 00242 " output 'bucket' is written to prefixname \n" 00243 "The mth sub-brick will contain: \n" 00244 "[-brick m coef k label] kth parameter regression coefficient \n" 00245 "[-brick m fstat label] F-stat for significance of regression \n" 00246 "[-brick m rstat label] coefficient of multiple determination R^2 \n" 00247 "[-brick m tstat k label] t-stat for kth regression coefficient \n" 00248 "\n" ); 00249 00250 printf 00251 ( 00252 "\n" 00253 "N.B.: For this program, the user must specify 1 and only 1 sub-brick \n" 00254 " with each -xydata command. That is, if an input dataset contains\n" 00255 " more than 1 sub-brick, a sub-brick selector must be used, e.g.: \n" 00256 " -xydata 2.17 4.59 7.18 'fred+orig[3]' \n" 00257 ); 00258 00259 printf("\n" MASTER_SHORTHELP_STRING ) ; 00260 00261 exit(0); 00262 } |
|
compute start time of this timeseries * Definition at line 2497 of file 3dRegAna.c. References EDIT_coerce_scale_type(), MCW_vol_amax(), MRI_IS_INT_TYPE, and top.
02505 { 02506 float fac=0.0 , top ; 02507 02508 if( MRI_IS_INT_TYPE(otype) ){ 02509 top = MCW_vol_amax( nxyz,1,1 , itype,ivol ) ; 02510 if (top == 0.0) fac = 0.0; 02511 else fac = MRI_TYPE_maxval[otype]/top; 02512 } 02513 02514 EDIT_coerce_scale_type( nxyz , fac , itype,ivol , otype,ovol ) ; 02515 return ( fac ); 02516 } |
|
Definition at line 312 of file 3dRegAna.c. References THD_3dim_dataset::daxes, ISVALID_3DIM_DATASET, THD_dataxes::nxx, THD_dataxes::nyy, THD_dataxes::nzz, THD_delete_3dim_dataset(), and THD_open_dataset().
00316 { 00317 00318 THD_3dim_dataset * dset=NULL; 00319 00320 /*----- read first dataset to get dimensions, etc. -----*/ 00321 00322 dset = THD_open_dataset( option_data->first_dataset ) ; 00323 if( ! ISVALID_3DIM_DATASET(dset) ){ 00324 fprintf(stderr,"*** Unable to open dataset file %s\n", 00325 option_data->first_dataset); 00326 exit(1) ; 00327 } 00328 00329 /*----- data set dimensions in voxels -----*/ 00330 option_data->nx = dset->daxes->nxx ; 00331 option_data->ny = dset->daxes->nyy ; 00332 option_data->nz = dset->daxes->nzz ; 00333 option_data->nxyz = option_data->nx * option_data->ny * option_data->nz ; 00334 00335 THD_delete_3dim_dataset( dset , False ) ; dset = NULL ; 00336 00337 } |
|
Definition at line 1077 of file 3dRegAna.c. References AFNI_logger(), argc, cols, display_help_menu(), FUNC_FIM_TYPE, FUNC_THR_TYPE, initialize_options(), ISVALID_3DIM_DATASET, malloc, matrix_create(), MAX_NAME_LENGTH, MAX_OBSERVATIONS, MAX_XVARS, p, PROGRAM_NAME, RA_error, rows, sort_model_indices(), THD_delete_3dim_dataset(), and THD_open_dataset(). Referenced by initialize_program().
01085 { 01086 const MAX_BRICKS = 100; /* max. number of bricks in the bucket */ 01087 int nopt = 1; /* input option argument counter */ 01088 int ival, index; /* integer input */ 01089 float fval; /* float input */ 01090 int rows, cols; /* number rows and columns for x matrix */ 01091 int irows, jcols; /* data point counters */ 01092 THD_3dim_dataset * dset=NULL; /* test whether data set exists */ 01093 char message[MAX_NAME_LENGTH]; /* error message */ 01094 char label[MAX_NAME_LENGTH]; /* sub-brick label */ 01095 int ibrick; /* sub-brick index */ 01096 int nbricks; /* number of sub-bricks in the bucket */ 01097 int p; /* number of parameters */ 01098 int ip, ix; /* parameter indices */ 01099 01100 01101 /*----- does user request help menu? -----*/ 01102 if (argc < 2 || strncmp(argv[1], "-help", 5) == 0) display_help_menu(); 01103 01104 01105 /*----- add to program log -----*/ 01106 AFNI_logger (PROGRAM_NAME,argc,argv); 01107 01108 01109 /*----- initialize the input options -----*/ 01110 initialize_options (regmodel, option_data); 01111 01112 01113 /*----- main loop over input options -----*/ 01114 while (nopt < argc && argv[nopt][0] == '-') 01115 { 01116 01117 /*----- -datum type -----*/ 01118 if( strncmp(argv[nopt],"-datum",6) == 0 ){ 01119 if( ++nopt >= argc ) RA_error("need an argument after -datum!") ; 01120 01121 if( strcmp(argv[nopt],"short") == 0 ){ 01122 option_data->datum = MRI_short ; 01123 } else if( strcmp(argv[nopt],"float") == 0 ){ 01124 option_data->datum = MRI_float ; 01125 } else { 01126 char buf[256] ; 01127 sprintf(buf,"-datum of type '%s' is not supported in 3dRegAna.", 01128 argv[nopt] ) ; 01129 RA_error(buf) ; 01130 } 01131 nopt++ ; continue ; /* go to next arg */ 01132 } 01133 01134 01135 /*----- -session dirname -----*/ 01136 if( strncmp(argv[nopt],"-session",8) == 0 ){ 01137 nopt++ ; 01138 if( nopt >= argc ) RA_error("need argument after -session!") ; 01139 strcpy(option_data->session , argv[nopt++]) ; 01140 continue ; 01141 } 01142 01143 01144 /*----- -diskspace -----*/ 01145 if( strncmp(argv[nopt],"-diskspace",10) == 0 ) 01146 { 01147 option_data->diskspace = 1; 01148 nopt++ ; continue ; /* go to next arg */ 01149 } 01150 01151 01152 /*----- -workmem megabytes -----*/ 01153 01154 if( strncmp(argv[nopt],"-workmem",6) == 0 ){ 01155 nopt++ ; 01156 if( nopt >= argc ) RA_error ("need argument after -workmem!") ; 01157 sscanf (argv[nopt], "%d", &ival); 01158 if( ival <= 0 ) RA_error ("illegal argument after -workmem!") ; 01159 option_data->workmem = ival ; 01160 nopt++ ; continue ; 01161 } 01162 01163 01164 /*----- -rmsmin r -----*/ 01165 if (strncmp(argv[nopt], "-rmsmin", 7) == 0) 01166 { 01167 nopt++; 01168 if (nopt >= argc) RA_error ("need argument after -rmsmin "); 01169 sscanf (argv[nopt], "%f", &fval); 01170 if (fval < 0.0) 01171 RA_error ("illegal argument after -rmsmin "); 01172 option_data->rms_min = fval; 01173 nopt++; 01174 continue; 01175 } 01176 01177 01178 /*----- -flof alpha -----*/ 01179 if (strncmp(argv[nopt], "-flof", 6) == 0) 01180 { 01181 nopt++; 01182 if (nopt >= argc) RA_error ("need argument after -flof "); 01183 sscanf (argv[nopt], "%f", &fval); 01184 if ((fval < 0.0) || (fval > 1.0)) 01185 RA_error ("illegal argument after -flof "); 01186 option_data->flofmax = fval; 01187 nopt++; 01188 continue; 01189 } 01190 01191 01192 /*----- -fdisp fval -----*/ 01193 if (strncmp(argv[nopt], "-fdisp", 6) == 0) 01194 { 01195 nopt++; 01196 if (nopt >= argc) RA_error ("need argument after -fdisp "); 01197 sscanf (argv[nopt], "%f", &fval); 01198 option_data->fdisp = fval; 01199 nopt++; 01200 continue; 01201 } 01202 01203 01204 /*----- -rows n -----*/ 01205 if (strncmp(argv[nopt], "-rows", 5) == 0) 01206 { 01207 nopt++; 01208 if (nopt >= argc) RA_error ("need argument after -rows "); 01209 sscanf (argv[nopt], "%d", &ival); 01210 if ((ival <= 0) || (ival > MAX_OBSERVATIONS)) 01211 RA_error ("illegal argument after -rows "); 01212 rows = ival; 01213 nopt++; 01214 continue; 01215 } 01216 01217 01218 /*----- -cols m -----*/ 01219 if (strncmp(argv[nopt], "-cols", 5) == 0) 01220 { 01221 nopt++; 01222 if (nopt >= argc) RA_error ("need argument after -cols "); 01223 sscanf (argv[nopt], "%d", &ival); 01224 if ((ival <= 0) || (ival > MAX_XVARS)) 01225 RA_error ("illegal argument after -cols "); 01226 cols = ival; 01227 nopt++; 01228 01229 matrix_create (rows, cols+1, xdata); 01230 irows = -1; 01231 continue; 01232 } 01233 01234 01235 /*----- -xydata x1 ... xm y -----*/ 01236 if (strncmp(argv[nopt], "-xydata", 7) == 0) 01237 { 01238 nopt++; 01239 if (nopt + cols >= argc) 01240 RA_error ("need cols+1 arguments after -xydata "); 01241 01242 irows++; 01243 if (irows > rows-1) RA_error ("too many data files"); 01244 01245 xdata->elts[irows][0] = 1.0; 01246 for (jcols = 1; jcols <= cols; jcols++) 01247 { 01248 sscanf (argv[nopt], "%f", &fval); 01249 xdata->elts[irows][jcols] = fval; 01250 nopt++; 01251 } 01252 01253 /*--- check whether input files exist ---*/ 01254 dset = THD_open_dataset( argv[nopt] ) ; 01255 if( ! ISVALID_3DIM_DATASET(dset) ) 01256 { 01257 sprintf(message,"Unable to open dataset file %s\n", argv[nopt]); 01258 RA_error (message); 01259 } 01260 THD_delete_3dim_dataset( dset , False ) ; dset = NULL ; 01261 01262 option_data->yname[irows] 01263 = malloc (sizeof(char) * MAX_NAME_LENGTH); 01264 strcpy (option_data->yname[irows], argv[nopt]); 01265 nopt++; 01266 continue; 01267 } 01268 01269 01270 /*----- -model -----*/ 01271 if (strncmp(argv[nopt], "-model", 6) == 0) 01272 { 01273 nopt++; 01274 if (nopt >= argc) RA_error ("need arguments after -model "); 01275 01276 while ((nopt < argc) 01277 && (strncmp(argv[nopt], "-", 1) != 0) 01278 && (strncmp(argv[nopt], ":", 1) != 0)) 01279 { 01280 sscanf (argv[nopt], "%d", &ival); 01281 if ((ival <= 0) || (ival > cols)) 01282 RA_error ("illegal argument after -model "); 01283 regmodel->flist[regmodel->p] = ival; 01284 01285 regmodel->p++; 01286 nopt++; 01287 } 01288 01289 if (strncmp(argv[nopt], ":", 1) == 0) 01290 { 01291 nopt++; 01292 01293 while ((nopt < argc) 01294 && (strncmp(argv[nopt], "-", 1) != 0)) 01295 { 01296 sscanf (argv[nopt], "%d", &ival); 01297 if ((ival < 0) || (ival > cols)) 01298 RA_error ("illegal argument after -model "); 01299 regmodel->flist[regmodel->p] = ival; 01300 regmodel->rlist[regmodel->q] = ival; 01301 regmodel->p++; 01302 regmodel->q++; 01303 nopt++; 01304 } 01305 } 01306 01307 /*----- sort model variable indices into ascending order -----*/ 01308 sort_model_indices (regmodel); 01309 01310 continue; 01311 } 01312 01313 01314 /*----- -fcoef k filename -----*/ 01315 if (strncmp(argv[nopt], "-fcoef", 6) == 0) 01316 { 01317 nopt++; 01318 if (nopt+1 >= argc) RA_error ("need 2 arguments after -fcoef "); 01319 sscanf (argv[nopt], "%d", &ival); 01320 if ((ival < 0) || (ival > cols)) 01321 RA_error ("illegal argument after -fcoef "); 01322 index = ival; 01323 nopt++; 01324 01325 option_data->fcoef_filename[index] = 01326 malloc (sizeof(char) * MAX_NAME_LENGTH); 01327 if (option_data->fcoef_filename[index] == NULL) 01328 RA_error ("Unable to allocate memory for fcoef_filename"); 01329 strcpy (option_data->fcoef_filename[index], argv[nopt]); 01330 01331 nopt++; 01332 continue; 01333 } 01334 01335 01336 /*----- -rcoef k filename -----*/ 01337 if (strncmp(argv[nopt], "-rcoef", 6) == 0) 01338 { 01339 nopt++; 01340 if (nopt+1 >= argc) RA_error ("need 2 arguments after -rcoef "); 01341 sscanf (argv[nopt], "%d", &ival); 01342 if ((ival < 0) || (ival > cols)) 01343 RA_error ("illegal argument after -rcoef "); 01344 index = ival; 01345 nopt++; 01346 01347 option_data->rcoef_filename[index] = 01348 malloc (sizeof(char) * MAX_NAME_LENGTH); 01349 if (option_data->rcoef_filename[index] == NULL) 01350 RA_error ("Unable to allocate memory for rcoef_filename"); 01351 strcpy (option_data->rcoef_filename[index], argv[nopt]); 01352 01353 nopt++; 01354 continue; 01355 } 01356 01357 01358 /*----- -tcoef k filename -----*/ 01359 if (strncmp(argv[nopt], "-tcoef", 6) == 0) 01360 { 01361 nopt++; 01362 if (nopt+1 >= argc) RA_error ("need 2 arguments after -tcoef "); 01363 sscanf (argv[nopt], "%d", &ival); 01364 if ((ival < 0) || (ival > cols)) 01365 RA_error ("illegal argument after -tcoef "); 01366 index = ival; 01367 nopt++; 01368 01369 option_data->tcoef_filename[index] = 01370 malloc (sizeof(char) * MAX_NAME_LENGTH); 01371 if (option_data->tcoef_filename[index] == NULL) 01372 RA_error ("Unable to allocate memory for tcoef_filename"); 01373 strcpy (option_data->tcoef_filename[index], argv[nopt]); 01374 01375 nopt++; 01376 continue; 01377 } 01378 01379 01380 /*----- -bucket n prefixname -----*/ 01381 if (strncmp(argv[nopt], "-bucket", 7) == 0) 01382 { 01383 nopt++; 01384 if (nopt+1 >= argc) RA_error ("need 2 arguments after -bucket "); 01385 sscanf (argv[nopt], "%d", &ival); 01386 if ((ival < 0) || (ival > MAX_BRICKS)) 01387 RA_error ("illegal argument after -bucket "); 01388 option_data->numbricks = ival; 01389 nopt++; 01390 01391 option_data->bucket_filename = 01392 malloc (sizeof(char) * MAX_NAME_LENGTH); 01393 if (option_data->bucket_filename == NULL) 01394 RA_error ("Unable to allocate memory for bucket_filename"); 01395 strcpy (option_data->bucket_filename, argv[nopt]); 01396 01397 /*----- set number of sub-bricks in the bucket -----*/ 01398 p = regmodel->p; 01399 if (ival == 0) 01400 nbricks = 2*p + 2; 01401 else 01402 nbricks = ival; 01403 option_data->numbricks = nbricks; 01404 01405 /*----- allocate memory and initialize bucket dataset options -----*/ 01406 option_data->brick_type = malloc (sizeof(int) * nbricks); 01407 option_data->brick_coef = malloc (sizeof(int) * nbricks); 01408 option_data->brick_label = malloc (sizeof(char *) * nbricks); 01409 for (ibrick = 0; ibrick < nbricks; ibrick++) 01410 { 01411 option_data->brick_type[ibrick] = -1; 01412 option_data->brick_coef[ibrick] = -1; 01413 option_data->brick_label[ibrick] = 01414 malloc (sizeof(char) * MAX_NAME_LENGTH); 01415 } 01416 01417 if (ival == 0) /*----- throw everything into the bucket -----*/ 01418 { 01419 for (ip = 0; ip < p; ip++) 01420 { 01421 ix = regmodel->flist[ip]; 01422 01423 ibrick = 2*ip; 01424 option_data->brick_type[ibrick] = FUNC_FIM_TYPE; 01425 option_data->brick_coef[ibrick] = ix; 01426 sprintf (label, "Coef #%d est.", ix); 01427 strcpy (option_data->brick_label[ibrick], label); 01428 01429 ibrick = 2*ip + 1; 01430 option_data->brick_type[ibrick] = FUNC_TT_TYPE; 01431 option_data->brick_coef[ibrick] = ix; 01432 sprintf (label, "Coef #%d t-stat", ix); 01433 strcpy (option_data->brick_label[ibrick], label); 01434 } 01435 01436 ibrick = 2*p; 01437 option_data->brick_type[ibrick] = FUNC_FT_TYPE; 01438 strcpy (option_data->brick_label[ibrick], "F-stat Regression"); 01439 01440 ibrick = 2*p + 1; 01441 option_data->brick_type[ibrick] = FUNC_THR_TYPE; 01442 strcpy (option_data->brick_label[ibrick], "R^2 Regression"); 01443 01444 } 01445 01446 nopt++; 01447 continue; 01448 } 01449 01450 01451 /*----- -brick m type k label -----*/ 01452 if (strncmp(argv[nopt], "-brick", 6) == 0) 01453 { 01454 nopt++; 01455 if (nopt+2 >= argc) RA_error ("need more arguments after -brick "); 01456 sscanf (argv[nopt], "%d", &ibrick); 01457 if ((ibrick < 0) || (ibrick >= option_data->numbricks)) 01458 RA_error ("illegal argument after -brick "); 01459 nopt++; 01460 01461 if (strncmp(argv[nopt], "coef", 4) == 0) 01462 { 01463 option_data->brick_type[ibrick] = FUNC_FIM_TYPE; 01464 01465 nopt++; 01466 sscanf (argv[nopt], "%d", &ival); 01467 if ((ival < 0) || (ival > cols)) 01468 RA_error ("illegal argument after coef "); 01469 option_data->brick_coef[ibrick] = ival; 01470 01471 nopt++; 01472 if (nopt >= argc) RA_error ("need more arguments after -brick "); 01473 strcpy (option_data->brick_label[ibrick], argv[nopt]); 01474 01475 } 01476 else if (strncmp(argv[nopt], "tstat", 4) == 0) 01477 { 01478 option_data->brick_type[ibrick] = FUNC_TT_TYPE; 01479 01480 nopt++; 01481 sscanf (argv[nopt], "%d", &ival); 01482 if ((ival < 0) || (ival > cols)) 01483 RA_error ("illegal argument after tstat "); 01484 option_data->brick_coef[ibrick] = ival; 01485 01486 nopt++; 01487 if (nopt >= argc) RA_error ("need more arguments after -brick "); 01488 strcpy (option_data->brick_label[ibrick], argv[nopt]); 01489 01490 } 01491 else if (strncmp(argv[nopt], "fstat", 4) == 0) 01492 { 01493 option_data->brick_type[ibrick] = FUNC_FT_TYPE; 01494 01495 nopt++; 01496 if (nopt >= argc) RA_error ("need more arguments after -brick "); 01497 strcpy (option_data->brick_label[ibrick], argv[nopt]); 01498 01499 } 01500 else if (strncmp(argv[nopt], "rstat", 4) == 0) 01501 { 01502 option_data->brick_type[ibrick] = FUNC_THR_TYPE; 01503 01504 nopt++; 01505 if (nopt >= argc) RA_error ("need more arguments after -brick "); 01506 strcpy (option_data->brick_label[ibrick], argv[nopt]); 01507 01508 } 01509 else 01510 { 01511 sprintf(message,"Unrecognized option after -brick: %s\n", 01512 argv[nopt]); 01513 RA_error (message); 01514 } 01515 01516 nopt++; 01517 continue; 01518 } 01519 01520 01521 /*----- unknown command -----*/ 01522 sprintf(message,"Unrecognized command line option: %s\n", argv[nopt]); 01523 RA_error (message); 01524 01525 } 01526 01527 /*----- check for fewer than expected datasets -----*/ 01528 if (irows < rows-1) 01529 RA_error ("Fewer than expected datasets were entered"); 01530 } |
|
Definition at line 1673 of file 3dRegAna.c. References cdff(), i, malloc, MTEST, p, q, and RA_error. Referenced by initialize_program().
01679 { 01680 int i, k; /* matrix row indices */ 01681 int j; /* matrix column index */ 01682 int match; /* boolean for repeat observation found */ 01683 01684 int which; /* indicator for F-distribution calculation */ 01685 double p, q; /* cumulative probabilities under F-dist. */ 01686 double f; /* F value corresponding to probability p */ 01687 double dfn, dfd; /* numerator and denominator dof */ 01688 int status; /* calculation status */ 01689 double bound; /* search bound */ 01690 01691 01692 /*----- allocate memory -----*/ 01693 option_data->levels = (int *) malloc (sizeof(int) * (xdata->rows)); 01694 MTEST (option_data->levels); 01695 option_data->counts = (int *) malloc (sizeof(int) * (xdata->rows)); 01696 MTEST (option_data->counts); 01697 01698 01699 /*----- initialization -----*/ 01700 for (i = 0; i < xdata->rows; i++) 01701 option_data->counts[i] = 0; 01702 option_data->levels[0] = 0; 01703 option_data->counts[0] = 1; 01704 option_data->c = 1; 01705 01706 01707 /*----- loop over observations -----*/ 01708 for (i = 1; i < xdata->rows; i++) 01709 { 01710 option_data->levels[i] = option_data->c; 01711 01712 /*----- determine if this is a repeat observation -----*/ 01713 for (k = 0; k < i; k++) 01714 { 01715 match = 1; 01716 for (j = 1; j < xdata->cols; j++) 01717 if (xdata->elts[i][j] != xdata->elts[k][j]) match = 0; 01718 01719 if (match) 01720 { 01721 option_data->levels[i] = option_data->levels[k]; 01722 break; 01723 } 01724 } 01725 01726 /*----- increment count of repeat observations -----*/ 01727 k = option_data->levels[i]; 01728 option_data->counts[k] ++; 01729 01730 /*----- increment count of distinct observation levels -----*/ 01731 if (k == option_data->c) option_data->c++; 01732 } 01733 01734 01735 /*----- determine F value corresponding to alpha for lack of fit -----*/ 01736 which = 2; 01737 q = option_data->flofmax; 01738 p = 1.0 - q; 01739 dfn = option_data->c - regmodel->p; 01740 dfd = xdata->rows - option_data->c; 01741 cdff (&which, &p, &q, &f, &dfn, &dfd, &status, &bound); 01742 if (status != 0) RA_error ("Error in calculating F - lack of fit "); 01743 option_data->flofmax = f; 01744 01745 01746 } |
|
Definition at line 2080 of file 3dRegAna.c. References calc_matrices(), matrix_destroy(), matrix_initialize(), MAX_XVARS, p, and q.
02087 : 1/(X'X) for full model */ 02088 matrix * xtxinvxt_full, /* matrix: (1/(X'X))X' for full model */ 02089 matrix * x_base, /* extracted X matrix for baseline model */ 02090 matrix * xtxinvxt_base, /* matrix: (1/(X'X))X' for baseline model */ 02091 matrix * x_rdcd, /* extracted X matrices for reduced models */ 02092 matrix * xtxinvxt_rdcd /* matrix: (1/(X'X))X' for reduced models */ 02093 ) 02094 02095 { 02096 int zlist[MAX_XVARS]; /* list of parameters in constant model */ 02097 int ip; /* parameter index */ 02098 matrix xtxinv_temp; /* intermediate results */ 02099 02100 02101 /*----- Initialize matrix -----*/ 02102 matrix_initialize (&xtxinv_temp); 02103 02104 02105 /*----- Initialize list of parameters in the constant model -----*/ 02106 for (ip = 0; ip < MAX_XVARS; ip++) 02107 zlist[ip] = 0; 02108 02109 02110 /*----- Calculate constant matrices which will be needed later -----*/ 02111 calc_matrices (xdata, 1, zlist, x_base, &xtxinv_temp, xtxinvxt_base); 02112 calc_matrices (xdata, q, rlist, x_rdcd, &xtxinv_temp, xtxinvxt_rdcd); 02113 calc_matrices (xdata, p, flist, x_full, xtxinv_full, xtxinvxt_full); 02114 02115 02116 /*----- Destroy matrix -----*/ 02117 matrix_destroy (&xtxinv_temp); 02118 02119 } |
|
Definition at line 943 of file 3dRegAna.c. References malloc, MAX_OBSERVATIONS, MAX_XVARS, and RA_error.
00948 { 00949 int ip; /* index */ 00950 00951 00952 regmodel->p = 0; 00953 regmodel->q = 0; 00954 00955 option_data->datum = ILLEGAL_TYPE; 00956 strcpy (option_data->session, "./"); 00957 00958 option_data->yname = NULL; 00959 option_data->first_dataset = NULL; 00960 00961 option_data->nx = 0; 00962 option_data->ny = 0; 00963 option_data->nz = 0; 00964 option_data->nxyz = 0; 00965 00966 option_data->diskspace = 0; 00967 option_data->workmem = 12; 00968 option_data->piece_size = 0; 00969 option_data->num_pieces = 0; 00970 00971 option_data->rms_min = 0.0; 00972 option_data->fdisp = -1.0; 00973 00974 option_data->levels = NULL; 00975 option_data->counts = NULL; 00976 option_data->c = 0; 00977 option_data->flofmax = -1; 00978 00979 option_data->numf = 0; 00980 option_data->numr = 0; 00981 option_data->numc = 0; 00982 option_data->numt = 0; 00983 00984 option_data->fcoef_filename = NULL; 00985 option_data->rcoef_filename = NULL; 00986 option_data->tcoef_filename = NULL; 00987 00988 option_data->numfiles = 0; 00989 00990 /*----- allocate memory for storing data file names -----*/ 00991 option_data->yname 00992 = (char **) malloc (sizeof(char *) * MAX_OBSERVATIONS); 00993 for (ip = 0; ip < MAX_OBSERVATIONS; ip++) 00994 option_data->yname[ip] = NULL; 00995 00996 00997 /*----- allocate memory space and initialize pointers for filenames -----*/ 00998 option_data->fcoef_filename = 00999 (char **) malloc (sizeof(char *) * MAX_XVARS); 01000 if (option_data->fcoef_filename == NULL) 01001 RA_error ("Unable to allocate memory for fcoef_filename"); 01002 01003 option_data->rcoef_filename = 01004 (char **) malloc (sizeof(char *) * MAX_XVARS); 01005 if (option_data->rcoef_filename == NULL) 01006 RA_error ("Unable to allocate memory for rcoef_filename"); 01007 01008 option_data->tcoef_filename = 01009 (char **) malloc (sizeof(char *) * MAX_XVARS); 01010 if (option_data->tcoef_filename == NULL) 01011 RA_error ("Unable to allocate memory for tcoef_filename"); 01012 01013 for (ip = 0; ip < MAX_XVARS; ip++) 01014 { 01015 option_data->fcoef_filename[ip] = NULL; 01016 option_data->rcoef_filename[ip] = NULL; 01017 option_data->tcoef_filename[ip] = NULL; 01018 } 01019 01020 01021 /*----- initialize bucket dataset options -----*/ 01022 option_data->bucket_filename = NULL; 01023 option_data->numbricks = -1; 01024 option_data->brick_type = NULL; 01025 option_data->brick_coef = NULL; 01026 option_data->brick_label = NULL; 01027 01028 } |
|
Definition at line 2014 of file 3dRegAna.c. References argc, break_into_pieces(), check_disk_space(), check_for_valid_inputs(), check_output_files(), check_temporary_files(), commandline, count_volumes_and_files(), get_dimensions(), get_inputs(), identify_repeats(), matrix_initialize(), PROGRAM_NAME, and tross_commandline().
02022 { 02023 02024 /*----- save command line for history notes -----*/ 02025 commandline = tross_commandline( PROGRAM_NAME , argc,argv ) ; 02026 02027 02028 /*----- create independent variable data matrix -----*/ 02029 matrix_initialize (xdata); 02030 02031 02032 /*----- get all operator inputs -----*/ 02033 get_inputs (argc, argv, xdata, regmodel, option_data); 02034 02035 02036 /*----- use first data set to get data set dimensions -----*/ 02037 option_data->first_dataset = option_data->yname[0]; 02038 get_dimensions (option_data); 02039 printf ("Data set dimensions: nx = %d ny = %d nz = %d nxyz = %d \n", 02040 option_data->nx, option_data->ny, option_data->nz, option_data->nxyz); 02041 02042 02043 /*----- count the number of output volumes and files required -----*/ 02044 count_volumes_and_files (regmodel, option_data); 02045 02046 02047 /*----- break problem into smaller pieces -----*/ 02048 break_into_pieces (xdata->rows, option_data); 02049 02050 02051 /*----- identify repeat observations -----*/ 02052 if (option_data->flofmax >= 0.0) 02053 identify_repeats (xdata, regmodel, option_data); 02054 02055 02056 /*----- check for valid inputs -----*/ 02057 check_for_valid_inputs (xdata, regmodel, option_data); 02058 02059 02060 /*----- check whether any of the output files already exist -----*/ 02061 check_output_files (option_data); 02062 02063 02064 /*----- check whether temporary files already exist -----*/ 02065 check_temporary_files (xdata, regmodel, option_data); 02066 02067 02068 /*----- check whether there is sufficient disk space -----*/ 02069 if (option_data->diskspace) check_disk_space (option_data); 02070 02071 } |
|
Definition at line 3184 of file 3dRegAna.c. References addto_args(), argc, calculate_results(), initialize_program(), machdep(), output_results(), PROGRAM_AUTHOR, PROGRAM_INITIAL, PROGRAM_LATEST, PROGRAM_NAME, and terminate_program().
03189 { 03190 RA_options option_data; /* user input options */ 03191 matrix xdata; /* independent variable matrix */ 03192 model regmodel; /* linear regression model */ 03193 int piece_size; /* number of voxels in dataset piece */ 03194 03195 03196 /*----- Identify software -----*/ 03197 printf ("\n\n"); 03198 printf ("Program: %s \n", PROGRAM_NAME); 03199 printf ("Author: %s \n", PROGRAM_AUTHOR); 03200 printf ("Initial Release: %s \n", PROGRAM_INITIAL); 03201 printf ("Latest Revision: %s \n", PROGRAM_LATEST); 03202 printf ("\n"); 03203 03204 /*-- 22 Feb 1999: addto the arglist, if user wants to --*/ 03205 03206 machdep() ; 03207 { int new_argc ; char ** new_argv ; 03208 addto_args( argc , argv , &new_argc , &new_argv ) ; 03209 if( new_argv != NULL ){ argc = new_argc ; argv = new_argv ; } 03210 } 03211 03212 /*----- program initialization -----*/ 03213 initialize_program (argc, argv, &xdata, ®model, &option_data); 03214 03215 03216 /*----- perform regression analysis -----*/ 03217 calculate_results (xdata, ®model, &option_data); 03218 03219 03220 /*----- write requested output files -----*/ 03221 output_results (xdata, ®model, &option_data); 03222 03223 03224 /*----- end of program -----*/ 03225 terminate_program (&xdata, ®model, &option_data); 03226 03227 exit(0); 03228 } |
|
Definition at line 2873 of file 3dRegAna.c. References free, FUNC_THR_TYPE, malloc, MAX_NAME_LENGTH, MTEST, p, q, read_volume(), write_afni_data(), and write_bucket_data().
02879 { 02880 int p; /* number of parameters in full model */ 02881 int q; /* number of parameters in reduced model */ 02882 int n; /* number of data points */ 02883 int nxyz; /* number of voxels */ 02884 int ip, ix; /* parameter index */ 02885 int numdof, dendof; /* numerator and denominator degrees of freedom */ 02886 int piece_size; /* number of voxels in dataset piece */ 02887 int num_pieces; /* dataset is divided into this many pieces */ 02888 float * volume1 = NULL; /* volume data for 1st sub-brick */ 02889 float * volume2 = NULL; /* volume data for 2nd sub-brick */ 02890 char filename[MAX_NAME_LENGTH]; /* name for temporary data file */ 02891 02892 02893 /*----- initialize local variables -----*/ 02894 p = regmodel->p; 02895 q = regmodel->q; 02896 n = xdata.rows; 02897 nxyz = option_data->nxyz; 02898 piece_size = option_data->piece_size; 02899 num_pieces = option_data->num_pieces; 02900 02901 02902 /*----- allocate memory space for volume data -----*/ 02903 volume1 = (float *) malloc (sizeof(float) * nxyz); 02904 MTEST (volume1); 02905 volume2 = (float *) malloc (sizeof(float) * nxyz); 02906 MTEST (volume2); 02907 02908 02909 /*----- write t-statistics data files -----*/ 02910 if (option_data->numt > 0) 02911 { 02912 numdof = n - p; 02913 dendof = 0; 02914 02915 for (ip = 0; ip < p; ip++) 02916 { 02917 ix = regmodel->flist[ip]; 02918 02919 if (option_data->tcoef_filename[ix] != NULL) 02920 { 02921 sprintf (filename, "coef.%d", ix); 02922 read_volume (filename, volume1, nxyz, piece_size, num_pieces); 02923 02924 sprintf (filename, "tcoef.%d", ix); 02925 read_volume (filename, volume2, nxyz, piece_size, num_pieces); 02926 02927 write_afni_data (option_data, 02928 option_data->tcoef_filename[ix], 02929 volume1, volume2, 02930 FUNC_TT_TYPE, numdof, dendof); 02931 } 02932 } 02933 } 02934 02935 02936 /*----- write R^2 data files -----*/ 02937 if (option_data->numr > 0) 02938 { 02939 strcpy (filename, "rsqr"); 02940 read_volume (filename, volume2, nxyz, piece_size, num_pieces); 02941 02942 for (ip = 0; ip < p; ip++) 02943 { 02944 ix = regmodel->flist[ip]; 02945 02946 if (option_data->rcoef_filename[ix] != NULL) 02947 { 02948 sprintf (filename, "coef.%d", ix); 02949 read_volume (filename, volume1, nxyz, piece_size, num_pieces); 02950 02951 write_afni_data (option_data, 02952 option_data->rcoef_filename[ix], 02953 volume1, volume2, 02954 FUNC_THR_TYPE, 0, 0); 02955 } 02956 } 02957 } 02958 02959 02960 /*----- write F-statistics data files -----*/ 02961 if (option_data->numf > 0) 02962 { 02963 strcpy (filename, "freg"); 02964 read_volume (filename, volume2, nxyz, piece_size, num_pieces); 02965 02966 numdof = p - q; 02967 dendof = n - p; 02968 02969 for (ip = 0; ip < p; ip++) 02970 { 02971 ix = regmodel->flist[ip]; 02972 02973 if (option_data->fcoef_filename[ix] != NULL) 02974 { 02975 sprintf (filename, "coef.%d", ix); 02976 read_volume (filename, volume1, nxyz, piece_size, num_pieces); 02977 02978 write_afni_data (option_data, 02979 option_data->fcoef_filename[ix], 02980 volume1, volume2, 02981 FUNC_FT_TYPE, numdof, dendof); 02982 } 02983 } 02984 } 02985 02986 02987 /*----- deallocate memory space for volume data -----*/ 02988 free (volume1); volume1 = NULL; 02989 free (volume2); volume2 = NULL; 02990 02991 02992 /*----- write the bucket dataset -----*/ 02993 if (option_data->numbricks > 0) 02994 write_bucket_data (xdata, regmodel, option_data); 02995 02996 02997 } |
|
Definition at line 182 of file 3dRegAna.c. References PROGRAM_NAME.
00183 { 00184 fprintf (stderr, "%s Error: %s \n", PROGRAM_NAME, message); 00185 exit(1); 00186 } |
|
Definition at line 347 of file 3dRegAna.c. References DOPEN, DSET_BRICK_FACTOR, DSET_BRICK_TYPE, DSET_PRINCIPAL_VALUE, EDIT_coerce_scale_type(), nz, SUB_POINTER, and THD_delete_3dim_dataset().
00355 { 00356 int iv; /* index number of intensity sub-brick */ 00357 THD_3dim_dataset * dset=NULL; /* data set pointer */ 00358 void * vfim = NULL; /* image data pointer */ 00359 int nx, ny, nz, nxyz; /* data set dimensions in voxels */ 00360 00361 nx = option_data->nx; 00362 ny = option_data->ny; 00363 nz = option_data->nz; 00364 nxyz = option_data->nxyz; 00365 00366 00367 /*----- read in the data -----*/ 00368 DOPEN (dset, filename) ; 00369 iv = DSET_PRINCIPAL_VALUE(dset) ; 00370 00371 /*----- convert it to floats (in ffim) -----*/ 00372 SUB_POINTER (dset, iv, fim_offset, vfim) ; 00373 EDIT_coerce_scale_type (piece_len, DSET_BRICK_FACTOR(dset,iv), 00374 DSET_BRICK_TYPE(dset,iv), vfim, /* input */ 00375 MRI_float ,ffim ) ; /* output */ 00376 00377 THD_delete_3dim_dataset( dset , False ) ; dset = NULL ; 00378 } |
|
Definition at line 416 of file 3dRegAna.c. References far, MAX_NAME_LENGTH, RA_error, and SUFFIX. Referenced by read_volume().
00422 { 00423 char sfilename[MAX_NAME_LENGTH]; /* piece file name */ 00424 char message[MAX_NAME_LENGTH]; /* error message */ 00425 FILE * far = NULL; /* floating point input file */ 00426 int count; /* number of data items read from disk */ 00427 00428 00429 /*----- input file name -----*/ 00430 strcpy (sfilename, filename); 00431 strcat (sfilename, SUFFIX); 00432 00433 /*----- open temporary data file for input -----*/ 00434 far = fopen (sfilename, "r"); 00435 if (far == NULL) 00436 { 00437 sprintf (message, "Unable to open temporary file %s", sfilename); 00438 RA_error (message); 00439 } 00440 00441 /*----- read 3d data file -----*/ 00442 count = fread (fin, sizeof(float), size, far); 00443 fclose (far); 00444 00445 /*----- error in reading file? -----*/ 00446 if (count != size) 00447 { 00448 sprintf (message, "Error in reading temporary file %s", sfilename); 00449 RA_error (message); 00450 } 00451 } |
|
Definition at line 869 of file 3dRegAna.c. References MAX_NAME_LENGTH, and read_piece(). Referenced by output_results(), and write_bucket_data().
00877 { 00878 int piece; /* piece index */ 00879 int piece_len; /* number of voxels in current piece */ 00880 int fim_offset; /* array offset to current piece */ 00881 char sfilename[MAX_NAME_LENGTH]; /* name for temporary data file */ 00882 00883 00884 /*----- loop over the temporary data file pieces -----*/ 00885 for (piece = 0; piece < num_pieces; piece++) 00886 { 00887 00888 /*----- current offset into volume -----*/ 00889 fim_offset = piece * piece_size; 00890 00891 /*----- size of current piece -----*/ 00892 if (piece < num_pieces-1) 00893 piece_len = piece_size; 00894 else 00895 piece_len = nxyz - fim_offset; 00896 00897 /*----- read in the piece data -----*/ 00898 sprintf (sfilename, "%s.p%d", filename, piece); 00899 read_piece (sfilename, volume + fim_offset, piece_len); 00900 00901 } /* loop over pieces */ 00902 00903 } |
|
Definition at line 2128 of file 3dRegAna.c. References c, calc_coef(), calc_flof(), calc_freg(), calc_rsqr(), calc_sse(), calc_sspe(), calc_tcoef(), p, q, vector_create(), vector_destroy(), and vector_initialize().
02133 : 1/(X'X) for full model */ 02134 matrix xtxinvxt_full, /* matrix: (1/(X'X))X' for full model */ 02135 matrix x_base, /* extracted X matrix for baseline model */ 02136 matrix xtxinvxt_base, /* matrix: (1/(X'X))X' for baseline model */ 02137 matrix x_rdcd, /* extracted X matrix for reduced model */ 02138 matrix xtxinvxt_rdcd, /* matrix: (1/(X'X))X' for reduced model */ 02139 vector y, /* vector of measured data */ 02140 float rms_min, /* minimum rms error to reject zero model */ 02141 int * levels, /* indices for repeat observations */ 02142 int * counts, /* number of observations at each level */ 02143 int c, /* number of unique rows in ind. var. matrix */ 02144 float flofmax, /* max. allowed F-stat due to lack of fit */ 02145 float * flof, /* F-statistic for lack of fit */ 02146 vector * coef_full, /* regression parameters */ 02147 vector * scoef_full, /* std. devs. for regression parameters */ 02148 vector * tcoef_full, /* t-statistics for regression parameters */ 02149 float * freg, /* F-statistic for the full regression model */ 02150 float * rsqr /* coeff. of multiple determination R^2 */ 02151 ) 02152 02153 { 02154 float sse_base; /* error sum of squares, baseline model */ 02155 float sse_rdcd; /* error sum of squares, reduced model */ 02156 float sse_full; /* error sum of squares, full model */ 02157 float sspe; /* pure error sum of squares */ 02158 vector coef_temp; /* intermediate results */ 02159 02160 02161 /*----- Initialization -----*/ 02162 vector_initialize (&coef_temp); 02163 02164 02165 /*----- Calculate regression coefficients for baseline model -----*/ 02166 calc_coef (xtxinvxt_base, y, &coef_temp); 02167 02168 02169 /*----- Calculate the error sum of squares for the baseline model -----*/ 02170 sse_base = calc_sse (x_base, coef_temp, y); 02171 02172 02173 /*----- Stop here if variation about baseline is sufficiently low -----*/ 02174 if (sqrt(sse_base/N) < rms_min) 02175 { 02176 vector_create (p, coef_full); 02177 vector_create (p, scoef_full); 02178 vector_create (p, tcoef_full); 02179 *freg = 0.0; 02180 *rsqr = 0.0; 02181 vector_destroy (&coef_temp); 02182 return; 02183 } 02184 02185 02186 /*----- Calculate regression coefficients for the full model -----*/ 02187 calc_coef (xtxinvxt_full, y, coef_full); 02188 02189 02190 /*----- Calculate the error sum of squares for the full model -----*/ 02191 sse_full = calc_sse (x_full, *coef_full, y); 02192 02193 02194 /*----- Calculate t-statistics for the regression coefficients -----*/ 02195 calc_tcoef (N, p, sse_full, xtxinv_full, 02196 *coef_full, scoef_full, tcoef_full); 02197 02198 02199 /*----- Test for lack of fit -----*/ 02200 if (flofmax > 0.0) 02201 { 02202 /*----- Calculate the pure error sum of squares -----*/ 02203 sspe = calc_sspe (y, levels, counts, c); 02204 02205 /*----- Calculate F-statistic for lack of fit -----*/ 02206 *flof = calc_flof (N, p, c, sse_full, sspe); 02207 02208 if (*flof > flofmax) 02209 { 02210 vector_create (p, coef_full); 02211 vector_create (p, scoef_full); 02212 vector_create (p, tcoef_full); 02213 *freg = 0.0; 02214 *rsqr = 0.0; 02215 vector_destroy (&coef_temp); 02216 return; 02217 } 02218 } 02219 else 02220 *flof = -1.0; 02221 02222 02223 /*----- Calculate regression coefficients for reduced model -----*/ 02224 calc_coef (xtxinvxt_rdcd, y, &coef_temp); 02225 02226 02227 /*----- Calculate the error sum of squares for the reduced model -----*/ 02228 sse_rdcd = calc_sse (x_rdcd, coef_temp, y); 02229 02230 02231 /*----- Calculate F-statistic for significance of the regression -----*/ 02232 *freg = calc_freg (N, p, q, sse_full, sse_rdcd); 02233 02234 02235 /*----- Calculate coefficient of multiple determination R^2 -----*/ 02236 *rsqr = calc_rsqr (sse_full, sse_base); 02237 02238 02239 /*----- Dispose of vector -----*/ 02240 vector_destroy (&coef_temp); 02241 02242 } |
|
Definition at line 704 of file 3dRegAna.c. References MAX_NAME_LENGTH, MAX_XVARS, and write_piece(). Referenced by calculate_results().
00714 { 00715 int ip; /* parameter index */ 00716 char filename[MAX_NAME_LENGTH]; /* name for temporary data file */ 00717 00718 00719 /*----- save piece containing F-statistics -----*/ 00720 if (freg_piece != NULL) 00721 { 00722 sprintf (filename, "freg.p%d", piece); 00723 write_piece (filename, freg_piece, piece_len); 00724 } 00725 00726 00727 /*----- save piece containing R^2 -----*/ 00728 if (rsqr_piece != NULL) 00729 { 00730 sprintf (filename, "rsqr.p%d", piece); 00731 write_piece (filename, rsqr_piece, piece_len); 00732 } 00733 00734 00735 /*----- save pieces containing regression coefficients -----*/ 00736 if (coef_piece != NULL) 00737 { 00738 for (ip = 0; ip < MAX_XVARS; ip++) 00739 if (coef_piece[ip] != NULL) 00740 { 00741 sprintf (filename, "coef.%d.p%d", ip, piece); 00742 write_piece (filename, coef_piece[ip], piece_len); 00743 } 00744 } 00745 00746 00747 /*----- save pieces containing t-statistics -----*/ 00748 if (tcoef_piece != NULL) 00749 { 00750 for (ip = 0; ip < MAX_XVARS; ip++) 00751 if (tcoef_piece[ip] != NULL) 00752 { 00753 sprintf (filename, "tcoef.%d.p%d", ip, piece); 00754 write_piece (filename, tcoef_piece[ip], piece_len); 00755 } 00756 } 00757 00758 } |
|
Definition at line 2251 of file 3dRegAna.c.
02268 { 02269 int ip; /* parameter index */ 02270 int ix; /* x-variable index */ 02271 02272 02273 /*----- save regression results into piece data -----*/ 02274 if (freg_piece != NULL) freg_piece[iv] = freg; 02275 if (rsqr_piece != NULL) rsqr_piece[iv] = rsqr; 02276 02277 02278 /*----- save regression parameter estimates -----*/ 02279 for (ip = 0; ip < regmodel->p; ip++) 02280 { 02281 ix = regmodel->flist[ip]; 02282 02283 if (coef_piece[ix] != NULL) coef_piece[ix][iv] = coef.elts[ip]; 02284 02285 if (tcoef_piece[ix] != NULL) tcoef_piece[ix][iv] = tcoef.elts[ip]; 02286 02287 } 02288 02289 02290 /*----- if so requested, display results for this voxel -----*/ 02291 if ((fdisp >= 0.0) && (freg >= fdisp)) 02292 { 02293 printf ("\n\nVoxel #%d: \n", iv); 02294 printf ("\nY data: \n"); 02295 for (ip = 0; ip < y.dim; ip++) 02296 printf ("Y[%d] = %f \n", ip, y.elts[ip]); 02297 02298 if (flof >= 0.0) printf ("\nF lack of fit = %f \n", flof); 02299 printf ("\nF regression = %f \n", freg); 02300 printf ("R-squared = %f \n", rsqr); 02301 02302 printf ("\nRegression Coefficients: \n"); 02303 for (ip = 0; ip < coef.dim; ip++) 02304 { 02305 ix = regmodel->flist[ip]; 02306 printf ("b[%d] = %f ", ix, coef.elts[ip]); 02307 printf ("t[%d] = %f \n", ix, tcoef.elts[ip]); 02308 } 02309 02310 } 02311 02312 } |
|
Definition at line 1037 of file 3dRegAna.c. Referenced by get_inputs().
01041 { 01042 int i, j, temp; /* model variable index numbers */ 01043 01044 01045 /*----- sort full model indices into ascending order -----*/ 01046 for (i = 0; i < regmodel->p - 1; i++) 01047 for (j = i+1; j < regmodel->p; j++) 01048 if (regmodel->flist[i] > regmodel->flist[j]) 01049 { 01050 temp = regmodel->flist[i]; 01051 regmodel->flist[i] = regmodel->flist[j]; 01052 regmodel->flist[j] = temp; 01053 } 01054 else if (regmodel->flist[i] == regmodel->flist[j]) 01055 RA_error ("Duplicate variable indices in model definition"); 01056 01057 01058 /*----- sort reduced model indices into ascending order -----*/ 01059 for (i = 0; i < regmodel->q - 1; i++) 01060 for (j = i+1; j < regmodel->q; j++) 01061 if (regmodel->rlist[i] > regmodel->rlist[j]) 01062 { 01063 temp = regmodel->rlist[i]; 01064 regmodel->rlist[i] = regmodel->rlist[j]; 01065 regmodel->rlist[j] = temp; 01066 } 01067 01068 } |
|
Definition at line 3006 of file 3dRegAna.c. References delete_volume(), free, matrix_destroy(), MAX_NAME_LENGTH, MAX_OBSERVATIONS, MAX_XVARS, and p.
03012 { 03013 int p; /* number of parameters in full model */ 03014 int ip, ix; /* parameter index */ 03015 int ib; /* sub-brick index */ 03016 int nxyz; /* number of voxels */ 03017 int piece_size; /* number of voxels in dataset piece */ 03018 int num_pieces; /* dataset is divided into this many pieces */ 03019 char filename[MAX_NAME_LENGTH]; /* name for temporary data file */ 03020 03021 03022 /*----- initialize local variables -----*/ 03023 p = regmodel->p; 03024 nxyz = option_data->nxyz; 03025 piece_size = option_data->piece_size; 03026 num_pieces = option_data->num_pieces; 03027 03028 03029 /*----- delete F-statistics data files -----*/ 03030 if (option_data->numf > 0) 03031 { 03032 strcpy (filename, "freg"); 03033 delete_volume (filename, nxyz, piece_size, num_pieces); 03034 } 03035 03036 03037 /*----- delete R^2 data files -----*/ 03038 if (option_data->numr > 0) 03039 { 03040 strcpy (filename, "rsqr"); 03041 delete_volume (filename, nxyz, piece_size, num_pieces); 03042 } 03043 03044 03045 /*----- delete t-statistics data files -----*/ 03046 if (option_data->numt > 0) 03047 { 03048 for (ip = 0; ip < p; ip++) 03049 { 03050 ix = regmodel->flist[ip]; 03051 sprintf (filename, "tcoef.%d", ix); 03052 delete_volume (filename, nxyz, piece_size, num_pieces); 03053 } 03054 } 03055 03056 03057 /*----- delete regression coefficients data files -----*/ 03058 if (option_data->numc > 0) 03059 { 03060 for (ip = 0; ip < p; ip++) 03061 { 03062 ix = regmodel->flist[ip]; 03063 sprintf (filename, "coef.%d", ix); 03064 delete_volume (filename, nxyz, piece_size, num_pieces); 03065 } 03066 } 03067 03068 03069 /*----- dispose of data matrix -----*/ 03070 matrix_destroy (xdata); 03071 03072 03073 /*----- deallocate memory -----*/ 03074 if (option_data->yname != NULL) 03075 { 03076 for (ip = 0; ip < MAX_OBSERVATIONS; ip++) 03077 { 03078 if (option_data->yname[ip] != NULL) 03079 { 03080 free (option_data->yname[ip]); 03081 option_data->yname[ip] = NULL; 03082 } 03083 } 03084 free (option_data->yname); 03085 option_data->yname = NULL; 03086 } 03087 03088 if (option_data->fcoef_filename != NULL) 03089 { 03090 for (ip = 0; ip < MAX_XVARS; ip++) 03091 { 03092 if (option_data->fcoef_filename[ip] != NULL) 03093 { 03094 free (option_data->fcoef_filename[ip]); 03095 option_data->fcoef_filename[ip] = NULL; 03096 } 03097 } 03098 free (option_data->fcoef_filename); 03099 option_data->fcoef_filename = NULL; 03100 } 03101 03102 if (option_data->rcoef_filename != NULL) 03103 { 03104 for (ip = 0; ip < MAX_XVARS; ip++) 03105 { 03106 if (option_data->rcoef_filename[ip] != NULL) 03107 { 03108 free (option_data->rcoef_filename[ip]); 03109 option_data->rcoef_filename[ip] = NULL; 03110 } 03111 } 03112 free (option_data->rcoef_filename); 03113 option_data->rcoef_filename = NULL; 03114 } 03115 03116 if (option_data->tcoef_filename != NULL) 03117 { 03118 for (ip = 0; ip < MAX_XVARS; ip++) 03119 { 03120 if (option_data->tcoef_filename[ip] != NULL) 03121 { 03122 free (option_data->tcoef_filename[ip]); 03123 option_data->tcoef_filename[ip] = NULL; 03124 } 03125 } 03126 free (option_data->tcoef_filename); 03127 option_data->tcoef_filename = NULL; 03128 } 03129 03130 if (option_data->levels != NULL) 03131 { 03132 free (option_data->levels); 03133 option_data->levels = NULL; 03134 } 03135 03136 if (option_data->counts != NULL) 03137 { 03138 free (option_data->counts); 03139 option_data->counts = NULL; 03140 } 03141 03142 03143 if (option_data->bucket_filename != NULL) 03144 { 03145 free (option_data->bucket_filename); 03146 option_data->bucket_filename = NULL; 03147 } 03148 03149 if (option_data->brick_type != NULL) 03150 { 03151 free (option_data->brick_type); 03152 option_data->brick_type = NULL; 03153 } 03154 03155 if (option_data->brick_coef != NULL) 03156 { 03157 free (option_data->brick_coef); 03158 option_data->brick_coef = NULL; 03159 } 03160 03161 if (option_data->brick_label != NULL) 03162 { 03163 for (ib = 0; ib < option_data->numbricks; ib++) 03164 { 03165 if (option_data->brick_label[ib] != NULL) 03166 { 03167 free (option_data->brick_label[ib]); 03168 option_data->brick_label[ib] = NULL; 03169 } 03170 } 03171 free (option_data->brick_label); 03172 option_data->brick_label = NULL; 03173 } 03174 03175 } |
|
Definition at line 2532 of file 3dRegAna.c. References ADN_brick_fac, ADN_datum_array, ADN_directory_name, ADN_func_type, ADN_label1, ADN_malloc_type, ADN_none, ADN_nvals, ADN_prefix, ADN_self_name, ADN_stat_aux, ADN_type, commandline, DATABLOCK_MEM_MALLOC, THD_3dim_dataset::dblk, THD_datablock::diskptr, DSET_BRICK, DSET_BRIKNAME, EDIT_coerce_autoscale_new(), EDIT_dset_items(), EDIT_empty_copy(), FUNC_FT_SCALE_SHORT, FUNC_THR_SCALE_SHORT, FUNC_THR_TYPE, FUNC_TT_SCALE_SHORT, GEN_FUNC_TYPE, HEAD_FUNC_TYPE, THD_diskptr::header_name, ISHEAD, ISVALID_3DIM_DATASET, malloc, MAX_STAT_AUX, mri_datum_size(), mri_fix_data_pointer(), RA_error, THD_delete_3dim_dataset(), THD_is_file(), THD_load_statistics(), THD_open_dataset(), THD_write_3dim_dataset(), top, tross_Append_History(), and tross_multi_Append_History().
02542 { 02543 int nxyz; /* number of voxels */ 02544 int ii; /* voxel index */ 02545 THD_3dim_dataset * dset=NULL; /* input afni data set pointer */ 02546 THD_3dim_dataset * new_dset=NULL; /* output afni data set pointer */ 02547 int ierror; /* number of errors in editing data */ 02548 int ibuf[32]; /* integer buffer */ 02549 float fbuf[MAX_STAT_AUX]; /* float buffer */ 02550 float fimfac; /* scale factor for short data */ 02551 int output_datum; /* data type for output data */ 02552 short * tsp = NULL; /* 2nd sub-brick data pointer */ 02553 void * vdif = NULL; /* 1st sub-brick data pointer */ 02554 float top, bot, func_scale_short; /* parameters for scaling data */ 02555 int top_ss, bot_ss; /* 2nd sub-brick value limits */ 02556 char label[80]; /* label for output file history */ 02557 02558 02559 /*----- initialize local variables -----*/ 02560 nxyz = option_data->nxyz; 02561 02562 /*----- read first dataset -----*/ 02563 dset = THD_open_dataset (option_data->first_dataset) ; 02564 if( ! ISVALID_3DIM_DATASET(dset) ){ 02565 fprintf(stderr,"*** Unable to open dataset file %s\n", 02566 option_data->first_dataset); exit(1) ; 02567 } 02568 02569 02570 /*-- make an empty copy of this dataset, for eventual output --*/ 02571 new_dset = EDIT_empty_copy( dset ) ; 02572 02573 02574 /*----- Record history of dataset -----*/ 02575 02576 sprintf (label, "Output prefix: %s", filename); 02577 if( commandline != NULL ) 02578 tross_multi_Append_History( new_dset , commandline,label,NULL ) ; 02579 else 02580 tross_Append_History ( new_dset, label); 02581 02582 02583 output_datum = MRI_short ; 02584 ibuf[0] = output_datum ; ibuf[1] = MRI_short ; 02585 02586 02587 ierror = EDIT_dset_items( new_dset , 02588 ADN_prefix , filename , 02589 ADN_label1 , filename , 02590 ADN_directory_name , option_data->session , 02591 ADN_self_name , filename , 02592 ADN_type , ISHEAD(dset) ? HEAD_FUNC_TYPE : 02593 GEN_FUNC_TYPE , 02594 ADN_func_type , func_type , 02595 ADN_nvals , FUNC_nvals[func_type] , 02596 ADN_datum_array , ibuf , 02597 ADN_malloc_type, DATABLOCK_MEM_MALLOC , 02598 ADN_none ) ; 02599 02600 if( ierror > 0 ){ 02601 fprintf(stderr, 02602 "*** %d errors in attempting to create output dataset!\n", ierror ) ; 02603 exit(1) ; 02604 } 02605 02606 if( THD_is_file(new_dset->dblk->diskptr->header_name) ){ 02607 fprintf(stderr, 02608 "*** Output dataset file %s already exists--cannot continue!\a\n", 02609 new_dset->dblk->diskptr->header_name ) ; 02610 exit(1) ; 02611 } 02612 02613 /*----- deleting exemplar dataset -----*/ 02614 THD_delete_3dim_dataset( dset , False ) ; dset = NULL ; 02615 02616 02617 /*----- allocate memory for output data -----*/ 02618 vdif = (void *) malloc( mri_datum_size(output_datum) * nxyz ) ; 02619 tsp = (short *) malloc( sizeof(short) * nxyz ) ; 02620 02621 /*----- attach bricks to new data set -----*/ 02622 mri_fix_data_pointer (vdif, DSET_BRICK(new_dset,0)); 02623 mri_fix_data_pointer (tsp, DSET_BRICK(new_dset,1)); 02624 02625 02626 /*----- convert data type to output specification -----*/ 02627 fimfac = EDIT_coerce_autoscale_new (nxyz, 02628 MRI_float, ffim, 02629 output_datum, vdif); 02630 if (fimfac != 0.0) fimfac = 1.0 / fimfac; 02631 02632 02633 top_ss = 32700; 02634 02635 if (func_type == FUNC_THR_TYPE) /* threshold */ 02636 { 02637 func_scale_short = FUNC_THR_SCALE_SHORT; 02638 bot_ss = 0; 02639 } 02640 else if (func_type == FUNC_TT_TYPE) /* t-statistic */ 02641 { 02642 func_scale_short = FUNC_TT_SCALE_SHORT; 02643 bot_ss = -top_ss; 02644 } 02645 else if (func_type == FUNC_FT_TYPE) /* F-statistic */ 02646 { 02647 func_scale_short = FUNC_FT_SCALE_SHORT; 02648 bot_ss = 0; 02649 } 02650 else 02651 RA_error ("Illegal ouput dataset function type"); 02652 02653 top = top_ss / func_scale_short; 02654 bot = bot_ss / func_scale_short; 02655 02656 02657 for (ii = 0; ii < nxyz; ii++) 02658 { 02659 if (ftr[ii] > top) 02660 tsp[ii] = top_ss; 02661 else if (ftr[ii] < bot) 02662 tsp[ii] = bot_ss; 02663 else 02664 tsp[ii] = (short) (func_scale_short * ftr[ii] + 0.5); 02665 } 02666 02667 02668 /*----- write afni data set -----*/ 02669 if (func_type == FUNC_THR_TYPE) /* threshold */ 02670 printf("----- Writing `fith' dataset "); 02671 else if (func_type == FUNC_TT_TYPE) /* t-statistic */ 02672 printf("----- Writing `fitt' dataset "); 02673 else if (func_type == FUNC_FT_TYPE) /* F-statistic */ 02674 printf("----- Writing `fift' dataset "); 02675 02676 printf("into %s\n", DSET_BRIKNAME(new_dset) ) ; 02677 02678 fbuf[0] = numdof; 02679 fbuf[1] = dendof; 02680 for( ii=2 ; ii < MAX_STAT_AUX ; ii++ ) fbuf[ii] = 0.0 ; 02681 (void) EDIT_dset_items( new_dset , ADN_stat_aux , fbuf , ADN_none ) ; 02682 02683 fbuf[0] = (output_datum == MRI_short && fimfac != 1.0 ) ? fimfac : 0.0 ; 02684 fbuf[1] = 1.0 / func_scale_short ; 02685 (void) EDIT_dset_items( new_dset , ADN_brick_fac , fbuf , ADN_none ) ; 02686 02687 THD_load_statistics( new_dset ) ; 02688 THD_write_3dim_dataset( NULL,NULL , new_dset , True ) ; 02689 02690 02691 /*----- deallocate memory -----*/ 02692 THD_delete_3dim_dataset( new_dset , False ) ; new_dset = NULL ; 02693 02694 } |
|
Definition at line 2703 of file 3dRegAna.c. References ADN_directory_name, ADN_func_type, ADN_malloc_type, ADN_none, ADN_ntt, ADN_nvals, ADN_prefix, ADN_type, commandline, DATABLOCK_MEM_MALLOC, delete_volume(), DSET_BRIKNAME, DSET_HEADNAME, EDIT_BRICK_FACTOR, EDIT_BRICK_LABEL, EDIT_BRICK_TO_FIFT, EDIT_BRICK_TO_FITT, EDIT_coerce_autoscale_new(), EDIT_dset_items(), EDIT_empty_copy(), EDIT_substitute_brick(), free, FUNC_BUCK_TYPE, FUNC_FIM_TYPE, FUNC_THR_TYPE, HEAD_FUNC_TYPE, malloc, MAX_NAME_LENGTH, MTEST, p, q, read_volume(), THD_delete_3dim_dataset(), THD_is_file(), THD_load_statistics(), THD_open_dataset(), THD_write_3dim_dataset(), and tross_Append_History().
02709 { 02710 const float EPSILON = 1.0e-10; 02711 int p; /* number of parameters in full model */ 02712 int q; /* number of parameters in reduced model */ 02713 int n; /* number of data points */ 02714 int nxyz; /* number of voxels */ 02715 THD_3dim_dataset * old_dset = NULL; /* prototype dataset */ 02716 THD_3dim_dataset * new_dset = NULL; /* output bucket dataset */ 02717 char * output_prefix; /* prefix name for bucket dataset */ 02718 char * output_session; /* directory for bucket dataset */ 02719 int nbricks, ib; /* number of sub-bricks in bucket dataset */ 02720 short ** bar = NULL; /* bar[ib] points to data for sub-brick #ib */ 02721 float factor; /* factor is new scale factor for sub-brick #ib */ 02722 int brick_type; /* indicates statistical type of sub-brick */ 02723 int brick_coef; /* regression coefficient index for sub-brick */ 02724 char * brick_label; /* character string label for sub-brick */ 02725 int ierror; /* number of errors in editing data */ 02726 char filename[MAX_NAME_LENGTH]; /* name for temporary data file */ 02727 int piece_size; /* number of voxels in dataset piece */ 02728 int num_pieces; /* dataset is divided into this many pieces */ 02729 float * volume = NULL; /* volume of floating point data */ 02730 char label[80]; /* label for output file history */ 02731 02732 02733 /*----- initialize local variables -----*/ 02734 p = regmodel->p; 02735 q = regmodel->q; 02736 n = xdata.rows; 02737 nxyz = option_data->nxyz; 02738 piece_size = option_data->piece_size; 02739 num_pieces = option_data->num_pieces; 02740 nbricks = option_data->numbricks; 02741 output_prefix = option_data->bucket_filename; 02742 output_session = option_data->session; 02743 02744 02745 /*----- allocate memory -----*/ 02746 volume = (float *) malloc (sizeof(float) * nxyz); 02747 MTEST (volume); 02748 bar = (short **) malloc (sizeof(short *) * nbricks); 02749 MTEST (bar); 02750 02751 02752 /*----- read first dataset -----*/ 02753 old_dset = THD_open_dataset (option_data->first_dataset) ; 02754 02755 02756 /*-- make an empty copy of this dataset, for eventual output --*/ 02757 new_dset = EDIT_empty_copy (old_dset); 02758 02759 02760 /*----- Record history of dataset -----*/ 02761 if( commandline != NULL ) 02762 tross_Append_History( new_dset , commandline ) ; 02763 sprintf (label, "Output prefix: %s", output_prefix); 02764 tross_Append_History ( new_dset, label); 02765 02766 02767 /*----- Modify some structural properties. Note that the nbricks 02768 just make empty sub-bricks, without any data attached. -----*/ 02769 ierror = EDIT_dset_items (new_dset, 02770 ADN_prefix, output_prefix, 02771 ADN_directory_name, output_session, 02772 ADN_type, HEAD_FUNC_TYPE, 02773 ADN_func_type, FUNC_BUCK_TYPE, 02774 ADN_ntt, 0, /* no time */ 02775 ADN_nvals, nbricks, 02776 ADN_malloc_type, DATABLOCK_MEM_MALLOC , 02777 ADN_none ) ; 02778 02779 if( ierror > 0 ) 02780 { 02781 fprintf(stderr, 02782 "*** %d errors in attempting to create output dataset!\n", 02783 ierror); 02784 exit(1); 02785 } 02786 02787 if (THD_is_file(DSET_HEADNAME(new_dset))) 02788 { 02789 fprintf(stderr, 02790 "*** Output dataset file %s already exists--cannot continue!\n", 02791 DSET_HEADNAME(new_dset)); 02792 exit(1); 02793 } 02794 02795 02796 /*----- deleting exemplar dataset -----*/ 02797 THD_delete_3dim_dataset( old_dset , False ); old_dset = NULL ; 02798 02799 02800 /*----- loop over new sub-brick index, attach data array with 02801 EDIT_substitute_brick then put some strings into the labels and 02802 keywords, and modify the sub-brick scaling factor -----*/ 02803 for (ib = 0; ib < nbricks; ib++) 02804 { 02805 /*----- get information about this sub-brick -----*/ 02806 brick_type = option_data->brick_type[ib]; 02807 brick_coef = option_data->brick_coef[ib]; 02808 brick_label = option_data->brick_label[ib]; 02809 02810 if (brick_type == FUNC_FIM_TYPE) 02811 { 02812 sprintf (filename, "coef.%d", brick_coef); 02813 } 02814 else if (brick_type == FUNC_THR_TYPE) 02815 { 02816 sprintf (filename, "rsqr"); 02817 } 02818 else if (brick_type == FUNC_TT_TYPE) 02819 { 02820 sprintf (filename, "tcoef.%d", brick_coef); 02821 EDIT_BRICK_TO_FITT (new_dset, ib, n-p); 02822 } 02823 else if (brick_type == FUNC_FT_TYPE) 02824 { 02825 sprintf (filename, "freg"); 02826 EDIT_BRICK_TO_FIFT (new_dset, ib, p-q, n-p); 02827 } 02828 02829 /*----- allocate memory for output sub-brick -----*/ 02830 bar[ib] = (short *) malloc (sizeof(short) * nxyz); 02831 MTEST (bar[ib]); 02832 02833 read_volume (filename, volume, nxyz, piece_size, num_pieces); 02834 delete_volume (filename, nxyz, piece_size, num_pieces); 02835 02836 factor = EDIT_coerce_autoscale_new (nxyz, MRI_float, volume, 02837 MRI_short, bar[ib]); 02838 02839 if (factor < EPSILON) factor = 0.0; 02840 else factor = 1.0 / factor; 02841 02842 /*----- edit the sub-brick -----*/ 02843 EDIT_BRICK_LABEL (new_dset, ib, brick_label); 02844 EDIT_BRICK_FACTOR (new_dset, ib, factor); 02845 02846 02847 /*----- attach bar[ib] to be sub-brick #ib -----*/ 02848 EDIT_substitute_brick (new_dset, ib, MRI_short, bar[ib]); 02849 02850 } 02851 02852 02853 /*----- write bucket data set -----*/ 02854 THD_load_statistics (new_dset); 02855 THD_write_3dim_dataset( NULL,NULL , new_dset , True ) ; 02856 fprintf(stderr,"----- Wrote bucket dataset "); 02857 fprintf(stderr,"into %s\n", DSET_BRIKNAME(new_dset)); 02858 02859 02860 /*----- deallocate memory -----*/ 02861 THD_delete_3dim_dataset( new_dset , False ) ; new_dset = NULL ; 02862 free (volume); volume = NULL; 02863 02864 } |
|
Definition at line 461 of file 3dRegAna.c. References far, fout, MAX_NAME_LENGTH, RA_error, and SUFFIX. Referenced by save_pieces().
00467 { 00468 char sfilename[MAX_NAME_LENGTH]; /* piece file name */ 00469 char message[MAX_NAME_LENGTH]; /* error message */ 00470 FILE * far = NULL; /* floating point output file */ 00471 int count; /* number of data items written to disk */ 00472 00473 00474 /*----- output file name -----*/ 00475 strcpy (sfilename, filename); 00476 strcat (sfilename, SUFFIX); 00477 00478 /*----- first, see if file already exists -----*/ 00479 far = fopen (sfilename, "r"); 00480 if (far != NULL) 00481 { 00482 sprintf (message, "Temporary file %s already exists. ", sfilename); 00483 RA_error (message); 00484 } 00485 00486 /*----- open temporary data file for output -----*/ 00487 far = fopen (sfilename, "w"); 00488 if (far == NULL) 00489 { 00490 sprintf (message, "Unable to write temporary file %s ", sfilename); 00491 RA_error (message); 00492 } 00493 00494 /*----- write 3d data set -----*/ 00495 count = fwrite (fout, sizeof(float), size, far); 00496 fclose (far); 00497 00498 /*----- error in writing file? -----*/ 00499 if (count != size) 00500 { 00501 sprintf (message, "Error in writing temporary file %s ", sfilename); 00502 RA_error (message); 00503 } 00504 00505 } |
Variable Documentation
|
Definition at line 120 of file 3dRegAna.c. Referenced by initialize_program(), write_afni_data(), and write_bucket_data(). |