Mercurial > hg > minc-tools
changeset 2384:d9c5d860f98c
* added znzlib and update nifti stuff. Now there is a seg fault that I need to
* fix
author | rotor <rotor> |
---|---|
date | Wed, 05 Sep 2007 06:22:36 +0000 (2007-09-05) |
parents | 1d2479d8da32 |
children | 16f2ab073f81 |
files | conversion/nifti1/mnc2nii.c conversion/nifti1/nifti1.h conversion/nifti1/nifti1_io.c conversion/nifti1/nifti1_io.h conversion/nifti1/nii2mnc.c conversion/nifti1/znzlib.c conversion/nifti1/znzlib.h |
diffstat | 7 files changed, 5473 insertions(+), 779 deletions(-) [+] |
line wrap: on
line diff
--- a/conversion/nifti1/mnc2nii.c +++ b/conversion/nifti1/mnc2nii.c @@ -109,7 +109,7 @@ nii_ptr->iname = NULL; nii_ptr->iname_offset = 0; nii_ptr->swapsize = 0; - nii_ptr->byteorder = LSB_FIRST; /* Use LSB first by default. */ + nii_ptr->byteorder = 1; /* default order (LSB) */ nii_ptr->data = NULL; for (i = 0; i < 4; i++) { @@ -556,7 +556,7 @@ nii_ptr->sto_xyz.m[3][2] = 0.0; nii_ptr->sto_xyz.m[3][3] = 1.0; - nii_ptr->sto_ijk = mat44_inverse(nii_ptr->sto_xyz); + nii_ptr->sto_ijk = nifti_mat44_inverse(nii_ptr->sto_xyz); nifti_datatype_sizes(nii_ptr->datatype, &nii_ptr->nbyper, &nii_ptr->swapsize); @@ -636,6 +636,7 @@ test_xform(nii_ptr->sto_xyz, 10, 10, 10); } + fprintf(stdout, "Calling NIFTI-1 Write routine\n"); nifti_image_write(nii_ptr); return (0);
--- a/conversion/nifti1/nifti1.h +++ b/conversion/nifti1/nifti1.h @@ -284,8 +284,9 @@ \brief Data structure defining the fields of a header extension. */ struct nifti1_extension { - int esize , ecode ; - char *edata ; + int esize ; /*!< size of extension, in bytes (must be multiple of 16) */ + int ecode ; /*!< extension code, one of the NIFTI_ECODE_ values */ + char * edata ; /*!< raw data, with no byte swapping */ } ; typedef struct nifti1_extension nifti1_extension ; @@ -785,11 +786,18 @@ #define NIFTI_INTENT_PVAL 22 - /*! Data is ln(p-value) (no params). */ + /*! Data is ln(p-value) (no params). + To be safe, a program should compute p = exp(-abs(this_value)). + The nifti_stats.c library returns this_value + as positive, so that this_value = -log(p). */ + #define NIFTI_INTENT_LOGPVAL 23 - /*! Data is log10(p-value) (no params). */ + /*! Data is log10(p-value) (no params). + To be safe, a program should compute p = pow(10.,-abs(this_value)). + The nifti_stats.c library returns this_value + as positive, so that this_value = -log10(p). */ #define NIFTI_INTENT_LOG10PVAL 24 @@ -1284,10 +1292,12 @@ if slice_duration is positive, indicates the timing pattern of the slice acquisition. The following codes are defined: - NIFTI_SLICE_SEQ_INC - NIFTI_SLICE_SEQ_DEC - NIFTI_SLICE_ALT_INC - NIFTI_SLICE_ALT_DEC + NIFTI_SLICE_SEQ_INC == sequential increasing + NIFTI_SLICE_SEQ_DEC == sequential decreasing + NIFTI_SLICE_ALT_INC == alternating increasing + NIFTI_SLICE_ALT_DEC == alternating decreasing + NIFTI_SLICE_ALT_INC2 == alternating increasing #2 + NIFTI_SLICE_ALT_DEC2 == alternating decreasing #2 { slice_start } = Indicates the start and end of the slice acquisition { slice_end } = pattern, when slice_code is nonzero. These values are present to allow for the possible addition of @@ -1297,22 +1307,34 @@ slice_end=dim[slice_dim]-1 are the correct values. For these values to be meaningful, slice_start must be non-negative and slice_end must be greater than - slice_start. + slice_start. Otherwise, they should be ignored. The following table indicates the slice timing pattern, relative to time=0 for the first slice acquired, for some sample cases. Here, dim[slice_dim]=7 (there are 7 slices, labeled 0..6), slice_duration=0.1, and slice_start=1, slice_end=5 (1 padded slice on each end). - slice - index SEQ_INC SEQ_DEC ALT_INC ALT_DEC - 6 -- n/a n/a n/a n/a n/a = not applicable - 5 -- 0.4 0.0 0.2 0.0 (slice time offset - 4 -- 0.3 0.1 0.4 0.3 doesn't apply to - 3 -- 0.2 0.2 0.1 0.1 slices outside range - 2 -- 0.1 0.3 0.3 0.4 slice_start..slice_end) - 1 -- 0.0 0.4 0.0 0.2 - 0 -- n/a n/a n/a n/a + slice + index SEQ_INC SEQ_DEC ALT_INC ALT_DEC ALT_INC2 ALT_DEC2 + 6 : n/a n/a n/a n/a n/a n/a n/a = not applicable + 5 : 0.4 0.0 0.2 0.0 0.4 0.2 (slice time offset + 4 : 0.3 0.1 0.4 0.3 0.1 0.0 doesn't apply to + 3 : 0.2 0.2 0.1 0.1 0.3 0.3 slices outside + 2 : 0.1 0.3 0.3 0.4 0.0 0.1 the range + 1 : 0.0 0.4 0.0 0.2 0.2 0.4 slice_start .. + 0 : n/a n/a n/a n/a n/a n/a slice_end) + + The SEQ slice_codes are sequential ordering (uncommon but not unknown), + either increasing in slice number or decreasing (INC or DEC), as + illustrated above. + + The ALT slice codes are alternating ordering. The 'standard' way for + these to operate (without the '2' on the end) is for the slice timing + to start at the edge of the slice_start .. slice_end group (at slice_start + for INC and at slice_end for DEC). For the 'ALT_*2' slice_codes, the + slice timing instead starts at the first slice in from the edge (at + slice_start+1 for INC2 and at slice_end-1 for DEC2). This latter + acquisition scheme is found on some Siemens scanners. The fields freq_dim, phase_dim, slice_dim are all squished into the single byte field dim_info (2 bits each, since the values for each field are @@ -1345,11 +1367,13 @@ of the slices @{ */ -#define NIFTI_SLICE_UNKNOWN 0 -#define NIFTI_SLICE_SEQ_INC 1 -#define NIFTI_SLICE_SEQ_DEC 2 -#define NIFTI_SLICE_ALT_INC 3 -#define NIFTI_SLICE_ALT_DEC 4 +#define NIFTI_SLICE_UNKNOWN 0 +#define NIFTI_SLICE_SEQ_INC 1 +#define NIFTI_SLICE_SEQ_DEC 2 +#define NIFTI_SLICE_ALT_INC 3 +#define NIFTI_SLICE_ALT_DEC 4 +#define NIFTI_SLICE_ALT_INC2 5 /* 05 May 2005: RWCox */ +#define NIFTI_SLICE_ALT_DEC2 6 /* 05 May 2005: RWCox */ /* @} */ /*---------------------------------------------------------------------------*/
--- a/conversion/nifti1/nifti1_io.c +++ b/conversion/nifti1/nifti1_io.c @@ -1,4 +1,6 @@ -#include "nifti1_io.h" /*** typedefs, prototypes, macros, etc. ***/ +#define _NIFTI1_IO_C_ + +#include "nifti1_io.h" /* typedefs, prototypes, macros, etc. */ /*****===================================================================*****/ /***** Sample functions to deal with NIFTI-1 and ANALYZE files *****/ @@ -14,11 +16,1002 @@ /***** incidental or otherwise, caused by any use of this document. *****/ /*****===================================================================*****/ +/** \file nifti1_io.c + \brief main collection of nifti1 i/o routines + - written by Bob Cox, SSCC NIMH + - revised by Mark Jenkinson, FMRIB + - revised by Rick Reynolds, SSCC, NIMH + - revised by Kate Fissell, University of Pittsburgh + + The library history can be viewed via "nifti_tool -nifti_hist". + <br>The library version can be viewed via "nifti_tool -nifti_ver". + */ + +/*! global history and version strings, for printing */ +static char * gni_history[] = +{ + "----------------------------------------------------------------------\n" + "history (of nifti library changes):\n" + "\n", + "0.0 August, 2003 [rwcox]\n" + " (Robert W Cox of the National Institutes of Health, SSCC/DIRP/NIMH)\n" + " - initial version\n" + "\n", + "0.1 July/August, 2004 [Mark Jenkinson]\n" + " (FMRIB Centre, University of Oxford, UK)\n" + " - Mainly adding low-level IO and changing things to allow gzipped\n" + " files to be read and written\n" + " - Full backwards compatability should have been maintained\n" + "\n", + "0.2 16 Nov 2004 [rickr]\n" + " (Rick Reynolds of the National Institutes of Health, SSCC/DIRP/NIMH)\n" + " - included Mark's changes in the AFNI distribution (including znzlib/)\n" + " (HAVE_ZLIB is commented out for the standard distribution)\n" + " - modified nifti_validfilename() and nifti_makebasename()\n" + " - added nifti_find_file_extension()\n" + "\n", + "0.3 3 Dec 2004 [rickr]\n" + " - note: header extensions are not yet checked for\n" + " - added formatted history as global string, for printing\n" + " - added nifti_disp_lib_hist(), to display the nifti library history\n" + " - added nifti_disp_lib_version(), to display the nifti library history\n", + " - re-wrote nifti_findhdrname()\n" + " o used nifti_find_file_extension()\n" + " o changed order of file tests (default is .nii, depends on input)\n" + " o free hdrname on failure\n" + " - made similar changes to nifti_findimgname()\n" + " - check for NULL return from nifti_findhdrname() calls\n", + " - removed most of ERREX() macros\n" + " - modified nifti_image_read()\n" + " o added debug info and error checking (on gni_debug > 0, only)\n" + " o fail if workingname is NULL\n" + " o check for failure to open header file\n" + " o free workingname on failure\n" + " o check for failure of nifti_image_load()\n" + " o check for failure of nifti_convert_nhdr2nim()\n", + " - changed nifti_image_load() to int, and check nifti_read_buffer return\n" + " - changed nifti_read_buffer() to fail on short read, and to count float\n" + " fixes (to print on debug)\n" + " - changed nifti_image_infodump to print to stderr\n" + " - updated function header comments, or moved comments above header\n" + " - removed const keyword\n" + " - added LNI_FERR() macro for error reporting on input files\n" + "\n", + "0.4 10 Dec 2004 [rickr] - added header extensions\n" + " - in nifti1_io.h:\n" + " o added num_ext and ext_list to the definition of nifti_image\n" + " o made many functions static (more to follow)\n" + " o added LNI_MAX_NIA_EXT_LEN, for max nifti_type 3 extension length\n", + " - added __DATE__ to version output in nifti_disp_lib_version()\n" + " - added nifti_disp_matrix_orient() to print orientation information\n" + " - added '.nia' as a valid file extension in nifti_find_file_extension()\n" + " - added much more debug output\n" + " - in nifti_image_read(), in the case of an ASCII header, check for\n" + " extensions after the end of the header\n", + " - added nifti_read_extensions() function\n" + " - added nifti_read_next_extension() function\n" + " - added nifti_add_exten_to_list() function\n" + " - added nifti_check_extension() function\n" + " - added nifti_write_extensions() function\n" + " - added nifti_extension_size() function\n" + " - in nifti_set_iname_offest():\n" + " o adjust offset by the extension size and the extender size\n", + " o fixed the 'ceiling modulo 16' computation\n" + " - in nifti_image_write_hdr_img2(): \n" + " o added extension writing\n" + " o check for NULL return from nifti_findimgname()\n" + " - include number of extensions in nifti_image_to_ascii() output\n" + " - in nifti_image_from_ascii():\n" + " o return bytes_read as a parameter, computed from the final spos\n" + " o extract num_ext from ASCII header\n" + "\n", + "0.5 14 Dec 2004 [rickr] - added sub-brick reading functions\n" + " - added nifti_brick_list type to nifti1_io.h, along with new prototypes\n" + " - added main nifti_image_read_bricks() function, with description\n" + " - added nifti_image_load_bricks() - library function (requires nim)\n" + " - added valid_nifti_brick_list() - library function\n" + " - added free_NBL() - library function\n", + " - added update_nifti_image_for_brick_list() for dimension update\n" + " - added nifti_load_NBL_bricks(), nifti_alloc_NBL_mem(),\n" + " nifti_copynsort() and force_positive() (static functions)\n" + " - in nifti_image_read(), check for failed load only if read_data is set\n" + " - broke most of nifti_image_load() into nifti_image_load_prep()\n" + "\n", + "0.6 15 Dec 2004 [rickr] - added sub-brick writing functionality\n" + " - in nifti1_io.h, removed znzlib directory from include - all nifti\n" + " library files are now under the nifti directory\n" + " - nifti_read_extensions(): print no offset warning for nifti_type 3\n" + " - nifti_write_all_data():\n" + " o pass nifti_brick_list * NBL, for optional writing\n" + " o if NBL, write each sub-brick, sequentially\n", + " - nifti_set_iname_offset(): case 1 must have sizeof() cast to int\n" + " - pass NBL to nifti_image_write_hdr_img2(), and allow NBL or data\n" + " - added nifti_image_write_bricks() wrapper for ...write_hdr_img2()\n" + " - included compression abilities\n" + "\n", + "0.7 16 Dec 2004 [rickr] - minor changes to extension reading\n" + "\n", + "0.8 21 Dec 2004 [rickr] - restrict extension reading, and minor changes\n" + " - in nifti_image_read(), compute bytes for extensions (see remaining)\n" + " - in nifti_read_extensions(), pass 'remain' as space for extensions,\n" + " pass it to nifti_read_next_ext(), and update for each one read \n" + " - in nifti_check_extension(), require (size <= remain)\n", + " - in update_nifti_image_brick_list(), update nvox\n" + " - in nifti_image_load_bricks(), make explicit check for nbricks <= 0\n" + " - in int_force_positive(), check for (!list)\n" + " - in swap_nifti_header(), swap sizeof_hdr, and reorder to struct order\n" + " - change get_filesize functions to signed ( < 0 is no file or error )\n", + " - in nifti_validfilename(), lose redundant (len < 0) check\n" + " - make print_hex_vals() static\n" + " - in disp_nifti_1_header, restrict string field widths\n" + "\n", + "0.9 23 Dec 2004 [rickr] - minor changes\n" + " - broke ASCII header reading out of nifti_image_read(), into new\n" + " functions has_ascii_header() and read_ascii_image()\n", + " - check image_read failure and znzseek failure\n" + " - altered some debug output\n" + " - nifti_write_all_data() now returns an int\n" + "\n", + "0.10 29 Dec 2004 [rickr]\n" + " - renamed nifti_valid_extension() to nifti_check_extension()\n" + " - added functions nifti_makehdrname() and nifti_makeimgname()\n" + " - added function valid_nifti_extensions()\n" + " - in nifti_write_extensions(), check for validity before writing\n", + " - rewrote nifti_image_write_hdr_img2():\n" + " o set write_data and leave_open flags from write_opts\n" + " o add debug print statements\n" + " o use nifti_write_ascii_image() for the ascii case\n" + " o rewrote the logic of all cases to be easier to follow\n", + " - broke out code as nifti_write_ascii_image() function\n" + " - added debug to top-level write functions, and free the znzFile\n" + " - removed unused internal function nifti_image_open()\n" + "\n", + "0.11 30 Dec 2004 [rickr] - small mods\n" + " - moved static function prototypes from header to C file\n" + " - free extensions in nifti_image_free()\n" + "\n", + "1.0 07 Jan 2005 [rickr] - INITIAL RELEASE VERSION\n" + " - added function nifti_set_filenames()\n" + " - added function nifti_read_header()\n" + " - added static function nhdr_looks_good()\n" + " - added static function need_nhdr_swap()\n" + " - exported nifti_add_exten_to_list symbol\n", + " - fixed #bytes written in nifti_write_extensions()\n" + " - only modify offset if it is too small (nifti_set_iname_offset)\n" + " - added nifti_type 3 to nifti_makehdrname and nifti_makeimgname\n" + " - added function nifti_set_filenames()\n" + "\n", + "1.1 07 Jan 2005 [rickr]\n" + " - in nifti_read_header(), swap if needed\n" + "\n", + "1.2 07 Feb 2005 [kate fissell c/o rickr] \n" + " - nifti1.h: added doxygen comments for main struct and #define groups\n" + " - nifti1_io.h: added doxygen comments for file and nifti_image struct\n" + " - nifti1_io.h: added doxygen comments for file and some functions\n" + " - nifti1_io.c: changed nifti_copy_nim_info to use memcpy\n" + "\n", + "1.3 09 Feb 2005 [rickr]\n" + " - nifti1.h: added doxygen comments for extension structs\n" + " - nifti1_io.h: put most #defines in #ifdef _NIFTI1_IO_C_ block\n" + " - added a doxygen-style description to every exported function\n" + " - added doxygen-style comments within some functions\n" + " - re-exported many znzFile functions that I had made static\n" + " - re-added nifti_image_open (sorry, Mark)\n" + " - every exported function now has 'nifti' in the name (19 functions)\n", + " - made sure every alloc() has a failure test\n" + " - added nifti_copy_extensions function, for use in nifti_copy_nim_info\n" + " - nifti_is_gzfile: added initial strlen test\n" + " - nifti_set_filenames: added set_byte_order parameter option\n" + " (it seems appropriate to set the BO when new files are associated)\n" + " - disp_nifti_1_header: prints to stdout (a.o.t. stderr), with fflush\n" + "\n", + "1.4 23 Feb 2005 [rickr] - sourceforge merge\n" + " - merged into the nifti_io CVS directory structure at sourceforge.net\n" + " - merged in 4 changes by Mark, and re-added his const keywords\n" + " - cast some pointers to (void *) for -pedantic compile option\n" + " - added nifti_free_extensions()\n" + "\n", + "1.5 02 Mar 2005 [rickr] - started nifti global options\n" + " - gni_debug is now g_opts.debug\n" + " - added validity check parameter to nifti_read_header\n" + " - need_nhdr_swap no longer does test swaps on the stack\n" + "\n", + "1.6 05 April 2005 [rickr] - validation and collapsed_image_read\n" + " - added nifti_read_collapsed_image(), an interface for reading partial\n" + " datasets, specifying a subset of array indices\n" + " - for read_collapsed_image, added static functions: rci_read_data(),\n" + " rci_alloc_mem(), and make_pivot_list()\n", + " - added nifti_nim_is_valid() to check for consistency (more to do)\n" + " - added nifti_nim_has_valid_dims() to do many dimensions tests\n" + "\n", + "1.7 08 April 2005 [rickr]\n" + " - added nifti_update_dims_from_array() - to update dimensions\n" + " - modified nifti_makehdrname() and nifti_makeimgname():\n" + " if prefix has a valid extension, use it (else make one up)\n" + " - added nifti_get_intlist - for making an array of ints\n" + " - fixed init of NBL->bsize in nifti_alloc_NBL_mem() {thanks, Bob}\n" + "\n", + "1.8 14 April 2005 [rickr]\n" + " - added nifti_set_type_from_names(), for nifti_set_filenames()\n" + " (only updates type if number of files does not match it)\n" + " - added is_valid_nifti_type(), just to be sure\n" + " - updated description of nifti_read_collapsed_image() for *data change\n" + " (if *data is already set, assume memory exists for results)\n" + " - modified rci_alloc_mem() to allocate only if *data is NULL\n" + "\n", + "1.9 19 April 2005 [rickr]\n" + " - added extension codes NIFTI_ECODE_COMMENT and NIFTI_ECODE_XCEDE\n" + " - added nifti_type codes NIFTI_MAX_ECODE and NIFTI_MAX_FTYPE\n" + " - added nifti_add_extension() {exported}\n" + " - added nifti_fill_extension() as a static function\n" + " - added nifti_is_valid_ecode() {exported}\n", + " - nifti_type values are now NIFTI_FTYPE_* file codes\n" + " - in nifti_read_extensions(), decrement 'remain' by extender size, 4\n" + " - in nifti_set_iname_offset(), case 1, update if offset differs\n" + " - only output '-d writing nifti file' if debug > 1\n" + "\n", + "1.10 10 May 2005 [rickr]\n" + " - files are read using ZLIB only if they end in '.gz'\n" + "\n", + "1.11 12 August 2005 [kate fissell]\n" + " - Kate's 0.2 release packaging, for sourceforge\n" + "\n", + "1.12 17 August 2005 [rickr] - comment (doxygen) updates\n" + " - updated comments for most functions (2 updates from Cinly Ooi)\n" + " - added nifti_type_and_names_match()\n" + "\n", + "1.12a 24 August 2005 [rickr] - remove all tabs from Clibs/*/*.[ch]\n", + "1.12b 25 August 2005 [rickr] - changes by Hans Johnson\n", + "1.13 25 August 2005 [rickr]\n", + " - finished changes by Hans for Insight\n" + " - added const in all appropraite parameter locations (30-40)\n" + " (any pointer referencing data that will not change)\n" + " - shortened all string constants below 509 character limit\n" + "1.14 28 October 2005 [HJohnson]\n", + " - use nifti_set_filenames() in nifti_convert_nhdr2nim()\n" + "1.15 02 November 2005 [rickr]\n", + " - added skip_blank_ext to nifti_global_options\n" + " - added nifti_set_skip_blank_ext(), to set option\n" + " - if skip_blank_ext and no extensions, do not read/write extender\n" + "1.16 18 November 2005 [rickr]\n", + " - removed any test or access of dim[i], i>dim[0]\n" + " - do not set pixdim for collapsed dims to 1.0, leave them as they are\n" + " - added magic and dim[i] tests in nifti_hdr_looks_good()\n" + " - added 2 size_t casts\n" + "1.17 22 November 2005 [rickr]\n", + " - in hdr->nim, for i > dim[0], pass 0 or 1, else set to 1\n" + "1.18 02 March 2006 [rickr]\n", + " - in nifti_alloc_NBL_mem(), fixed nt=0 case from 1.17 change\n" + "1.19 23 May 2006 [HJohnson,rickr]\n", + " - nifti_write_ascii_image(): free(hstr)\n" + " - nifti_copy_extensions(): clear num_ext and ext_list\n" + "1.20 27 Jun 2006 [rickr]\n", + " - nifti_findhdrname(): fixed assign of efirst to match stated logic\n" + " (problem found by Atle Bjørnerud)\n" + "1.21 5 Sep 2006 [rickr] update for nifticlib-0.4 release\n", + " - was reminded to actually add nifti_set_skip_blank_ext()\n" + " - init g_opts.skip_blank_ext to 0\n" + "----------------------------------------------------------------------\n" +}; +static char gni_version[] = "nifti library version 1.21 5 Sep, 2006)"; + +/*! global nifti options structure */ +static nifti_global_options g_opts = { 1, 0 }; + /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI datatype. - Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/* prototypes for internal functions - not part of exported library */ + +/* extension routines */ +static int nifti_read_extensions( nifti_image *nim, znzFile fp, int remain ); +static int nifti_read_next_extension( nifti1_extension * nex, nifti_image *nim, int remain, znzFile fp ); +static int nifti_check_extension(nifti_image *nim, int size,int code, int rem); +static void update_nifti_image_for_brick_list(nifti_image * nim , int nbricks); +static int nifti_add_exten_to_list(nifti1_extension * new_ext, + nifti1_extension ** list, int new_length); +static int nifti_fill_extension(nifti1_extension * ext, const char * data, + int len, int ecode); + +/* NBL routines */ +static int nifti_load_NBL_bricks(nifti_image * nim , int * slist, int * sindex, nifti_brick_list * NBL, znzFile fp ); +static int nifti_alloc_NBL_mem( nifti_image * nim, int nbricks, + nifti_brick_list * nbl); +static int nifti_copynsort(int nbricks, const int *blist, int **slist, + int **sindex); + +/* for nifti_read_collapsed_image: */ +static int rci_read_data(nifti_image *nim, int *pivots, int *prods, int nprods, + const int dims[], char *data, znzFile fp, int base_offset); +static int rci_alloc_mem(void ** data, int prods[8], int nprods, int nbyper ); +static int make_pivot_list(nifti_image * nim, const int dims[], int pivots[], + int prods[], int * nprods ); + +/* misc */ +static int need_nhdr_swap (short dim0, int hdrsize); +static int print_hex_vals (const char * data, int nbytes, FILE * fp); +static int unescape_string (char *str); /* string utility functions */ +static char *escapize_string (const char *str); + +/* internal I/O routines */ +static znzFile nifti_image_load_prep( nifti_image *nim ); +static int has_ascii_header(znzFile fp); +/*---------------------------------------------------------------------------*/ + + +/* for calling from some main program */ + +/*----------------------------------------------------------------------*/ +/*! display the nifti library module history (via stdout) +*//*--------------------------------------------------------------------*/ +void nifti_disp_lib_hist( void ) +{ + int c, len = sizeof(gni_history)/sizeof(char *); + for( c = 0; c < len; c++ ) + fputs(gni_history[c], stdout); +} + +/*----------------------------------------------------------------------*/ +/*! display the nifti library version (via stdout) +*//*--------------------------------------------------------------------*/ +void nifti_disp_lib_version( void ) +{ + printf("%s, compiled %s\n", gni_version, __DATE__); +} + + +/*----------------------------------------------------------------------*/ +/*! nifti_image_read_bricks - read nifti data as array of bricks + * + * 13 Dec 2004 [rickr] + * + * \param hname - filename of dataset to read (must be valid) + * \param nbricks - number of sub-bricks to read + * (if blist is valid, nbricks must be > 0) + * \param blist - list of sub-bricks to read + * (can be NULL; if NULL, read complete dataset) + * \param NBL - pointer to empty nifti_brick_list struct + * (must be a valid pointer) + * + * \return + * <br> nim - same as nifti_image_read, but nim->data will be NULL + * <br> NBL - filled with data + * + * By default, this function will read the nifti dataset and break the data + * into a list of nt*nu*nv*nw sub-bricks, each having size nx*ny*nz elements. + * That is to say, instead of reading the entire dataset as a single array, + * break it up into sub-bricks, each of size nx*ny*nz elements. + * + * If 'blist' is valid, it is taken to be a list of sub-bricks, of length + * 'nbricks'. The data will still be separated into sub-bricks of size + * nx*ny*nz elements, but now 'nbricks' sub-bricks will be returned, of the + * caller's choosing via 'blist'. + * + * E.g. consider a dataset with 12 sub-bricks (numbered 0..11), and the + * following code: + * + * <pre> + * { nifti_brick_list NB_orig, NB_select; + * nifti_image * nim_orig, * nim_select; + * int blist[5] = { 7, 0, 5, 5, 9 }; + * + * nim_orig = nifti_image_read_bricks("myfile.nii", 0, NULL, &NB_orig); + * nim_select = nifti_image_read_bricks("myfile.nii", 5, blist, &NB_select); + * } + * </pre> + * + * Here, nim_orig gets the entire dataset, where NB_orig.nbricks = 11. But + * nim_select has NB_select.nbricks = 5. + * + * Note that the first case is not quite the same as just calling the + * nifti_image_read function, as here the data is separated into sub-bricks. + * + * Note that valid blist elements are in [0..nt*nu*nv*nw-1], + * or written [ 0 .. (dim[4]*dim[5]*dim[6]*dim[7] - 1) ]. + * + * Note that, as is the case with all of the reading functions, the + * data will be allocated, read in, and properly byte-swapped, if + * necessary. + * + * \sa nifti_image_load_bricks, nifti_free_NBL, valid_nifti_brick_list, + nifti_image_read +*//*----------------------------------------------------------------------*/ +nifti_image *nifti_image_read_bricks(const char * hname, int nbricks, + const int * blist, nifti_brick_list * NBL) +{ + nifti_image * nim; + + if( !hname || !NBL ){ + fprintf(stderr,"** nifti_image_read_bricks: bad params (%p,%p)\n", + hname, (void *)NBL); + return NULL; + } + + if( blist && nbricks <= 0 ){ + fprintf(stderr,"** nifti_image_read_bricks: bad nbricks, %d\n", nbricks); + return NULL; + } + + nim = nifti_image_read(hname, 0); /* read header, but not data */ + + if( !nim ) return NULL; /* errors were already printed */ + + /* if we fail, free image and return */ + if( nifti_image_load_bricks(nim, nbricks, blist, NBL) <= 0 ){ + nifti_image_free(nim); + return NULL; + } + + if( blist ) update_nifti_image_for_brick_list(nim, nbricks); + + return nim; +} + + +/*---------------------------------------------------------------------- + * update_nifti_image_for_brick_list - update nifti_image + * + * When loading a specific brick list, the distinction between + * nt, nu, nv and nw is lost. So put everything in t, and set + * dim[0] = 4. + *----------------------------------------------------------------------*/ +static void update_nifti_image_for_brick_list( nifti_image * nim , int nbricks ) +{ + int ndim; + + if( g_opts.debug > 2 ){ + fprintf(stderr,"+d updating image dimensions for %d bricks in list\n", + nbricks); + fprintf(stderr," ndim = %d\n",nim->ndim); + fprintf(stderr," nx,ny,nz,nt,nu,nv,nw: (%d,%d,%d,%d,%d,%d,%d)\n", + nim->nx, nim->ny, nim->nz, nim->nt, nim->nu, nim->nv, nim->nw); + } + + nim->nt = nbricks; + nim->nu = nim->nv = nim->nw = 1; + nim->dim[4] = nbricks; + nim->dim[5] = nim->dim[6] = nim->dim[7] = 1; + + /* compute nvox */ + /* do not rely on dimensions above dim[0] 16 Nov 2005 [rickr] */ + for( nim->nvox = 1, ndim = 1; ndim <= nim->dim[0]; ndim++ ) + nim->nvox *= nim->dim[ndim]; + + /* update the dimensions to 4 or lower */ + for( ndim = 4; (ndim > 1) && (nim->dim[ndim] <= 1); ndim-- ) + ; + + if( g_opts.debug > 2 ){ + fprintf(stderr,"+d ndim = %d -> %d\n",nim->ndim, ndim); + fprintf(stderr," --> (%d,%d,%d,%d,%d,%d,%d)\n", + nim->nx, nim->ny, nim->nz, nim->nt, nim->nu, nim->nv, nim->nw); + } + + nim->dim[0] = nim->ndim = ndim; +} + + +/*----------------------------------------------------------------------*/ +/*! nifti_update_dims_from_array - update nx, ny, ... from nim->dim[] + + Fix all the dimension information, based on a new nim->dim[]. + + Note: we assume that dim[0] will not increase. + + Check for updates to pixdim[], dx,..., nx,..., nvox, ndim, dim[0]. +*//*--------------------------------------------------------------------*/ +int nifti_update_dims_from_array( nifti_image * nim ) +{ + int c, ndim; + + if( !nim ){ + fprintf(stderr,"** update_dims: missing nim\n"); + return 1; + } + + if( g_opts.debug > 2 ){ + fprintf(stderr,"+d updating image dimensions given nim->dim:"); + for( c = 0; c < 8; c++ ) fprintf(stderr," %d", nim->dim[c]); + fputc('\n',stderr); + } + + /* verify dim[0] first */ + if(nim->dim[0] < 1 || nim->dim[0] > 7){ + fprintf(stderr,"** invalid dim[0], dim[] = "); + for( c = 0; c < 8; c++ ) fprintf(stderr," %d", nim->dim[c]); + fputc('\n',stderr); + return 1; + } + + /* set nx, ny ..., dx, dy, ..., one by one */ + + /* less than 1, set to 1, else copy */ + if(nim->dim[1] < 1) nim->nx = nim->dim[1] = 1; + else nim->nx = nim->dim[1]; + nim->dx = nim->pixdim[1]; + + /* if undefined, or less than 1, set to 1 */ + if(nim->dim[0] < 2 || (nim->dim[0] >= 2 && nim->dim[2] < 1)) + nim->ny = nim->dim[2] = 1; + else + nim->ny = nim->dim[2]; + /* copy delta values, in any case */ + nim->dy = nim->pixdim[2]; + + if(nim->dim[0] < 3 || (nim->dim[0] >= 3 && nim->dim[3] < 1)) + nim->nz = nim->dim[3] = 1; + else /* just copy vals from arrays */ + nim->nz = nim->dim[3]; + nim->dz = nim->pixdim[3]; + + if(nim->dim[0] < 4 || (nim->dim[0] >= 4 && nim->dim[4] < 1)) + nim->nt = nim->dim[4] = 1; + else /* just copy vals from arrays */ + nim->nt = nim->dim[4]; + nim->dt = nim->pixdim[4]; + + if(nim->dim[0] < 5 || (nim->dim[0] >= 5 && nim->dim[5] < 1)) + nim->nu = nim->dim[5] = 1; + else /* just copy vals from arrays */ + nim->nu = nim->dim[5]; + nim->du = nim->pixdim[5]; + + if(nim->dim[0] < 6 || (nim->dim[0] >= 6 && nim->dim[6] < 1)) + nim->nv = nim->dim[6] = 1; + else /* just copy vals from arrays */ + nim->nv = nim->dim[6]; + nim->dv = nim->pixdim[6]; + + if(nim->dim[0] < 7 || (nim->dim[0] >= 7 && nim->dim[7] < 1)) + nim->nw = nim->dim[7] = 1; + else /* just copy vals from arrays */ + nim->nw = nim->dim[7]; + nim->dw = nim->pixdim[7]; + + for( c = 1, nim->nvox = 1; c <= nim->dim[0]; c++ ) + nim->nvox *= nim->dim[c]; + + /* compute ndim, assuming it can be no larger than the old one */ + for( ndim = nim->dim[0]; (ndim > 1) && (nim->dim[ndim] <= 1); ndim-- ) + ; + + if( g_opts.debug > 2 ){ + fprintf(stderr,"+d ndim = %d -> %d\n",nim->ndim, ndim); + fprintf(stderr," --> (%d,%d,%d,%d,%d,%d,%d)\n", + nim->nx, nim->ny, nim->nz, nim->nt, nim->nu, nim->nv, nim->nw); + } + + nim->dim[0] = nim->ndim = ndim; + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! Load the image data from disk into an already-prepared image struct. + * + * \param nim - initialized nifti_image, without data + * \param nbricks - the length of blist (must be 0 if blist is NULL) + * \param blist - an array of xyz volume indices to read (can be NULL) + * \param NBL - pointer to struct where resulting data will be stored + * + * If blist is NULL, read all sub-bricks. + * + * \return the number of loaded bricks (NBL->nbricks), + * 0 on failure, < 0 on error + * + * NOTE: it is likely that another function will copy the data pointers + * out of NBL, in which case the only pointer the calling function + * will want to free is NBL->bricks (not each NBL->bricks[i]). +*//*--------------------------------------------------------------------*/ +int nifti_image_load_bricks( nifti_image * nim , int nbricks, + const int * blist, nifti_brick_list * NBL ) +{ + int * slist = NULL, * sindex = NULL, rv; + znzFile fp; + + /* we can have blist == NULL */ + if( !nim || !NBL ){ + fprintf(stderr,"** nifti_image_load_bricks, bad params (%p,%p)\n", + (void *)nim, (void *)NBL); + return -1; + } + + if( blist && nbricks <= 0 ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d load_bricks: received blist with nbricks = %d," + "ignoring blist\n", nbricks); + blist = NULL; /* pretend nothing was passed */ + } + + if( blist && ! valid_nifti_brick_list( nim, nbricks, blist, g_opts.debug>0 ) ) + return -1; + + /* for efficiency, let's read the file in order */ + if( blist && nifti_copynsort( nbricks, blist, &slist, &sindex ) != 0 ) + return -1; + + /* open the file and position the FILE pointer */ + fp = nifti_image_load_prep( nim ); + if( !fp ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** nifti_image_load_bricks, failed load_prep\n"); + if( blist ){ free(slist); free(sindex); } + return -1; + } + + /* this will flag to allocate defaults */ + if( !blist ) nbricks = 0; + if( nifti_alloc_NBL_mem( nim, nbricks, NBL ) != 0 ){ + if( blist ){ free(slist); free(sindex); } + znzclose(fp); + return -1; + } + + rv = nifti_load_NBL_bricks(nim, slist, sindex, NBL, fp); + + if( rv != 0 ){ + nifti_free_NBL( NBL ); /* failure! */ + NBL->nbricks = 0; /* repetative, but clear */ + } + + if( slist ){ free(slist); free(sindex); } + + znzclose(fp); + + return NBL->nbricks; +} + + +/*----------------------------------------------------------------------*/ +/*! nifti_free_NBL - free all pointers and clear structure + * + * note: this does not presume to free the structure pointer +*//*--------------------------------------------------------------------*/ +void nifti_free_NBL( nifti_brick_list * NBL ) +{ + int c; + + if( NBL->bricks ){ + for( c = 0; c < NBL->nbricks; c++ ) + if( NBL->bricks[c] ) free(NBL->bricks[c]); + free(NBL->bricks); + NBL->bricks = NULL; + } + + NBL->nbricks = NBL->bsize = 0; +} + + +/*---------------------------------------------------------------------- + * nifti_load_NBL_bricks - read the file data into the NBL struct + * + * return 0 on success, -1 on failure + *----------------------------------------------------------------------*/ +static int nifti_load_NBL_bricks( nifti_image * nim , int * slist, int * sindex, + nifti_brick_list * NBL, znzFile fp ) +{ + int c, rv; + int oposn, fposn; /* orig and current file positions */ + int prev, isrc, idest; /* previous and current sub-brick, and new index */ + + oposn = znztell(fp); /* store current file position */ + fposn = oposn; + if( fposn < 0 ){ + fprintf(stderr,"** load bricks: ztell failed??\n"); + return -1; + } + + /* first, handle the default case, no passed blist */ + if( !slist ){ + for( c = 0; c < NBL->nbricks; c++ ) { + rv = nifti_read_buffer(fp, NBL->bricks[c], NBL->bsize, nim); + if( rv != NBL->bsize ){ + fprintf(stderr,"** load bricks: cannot read brick %d from '%s'\n", + c, nim->iname ? nim->iname : nim->fname); + return -1; + } + } + if( g_opts.debug > 1 ) + fprintf(stderr,"+d read %d default bricks from file %s\n", + NBL->nbricks, nim->iname ? nim->iname : nim->fname ); + return 0; + } + + if( !sindex ){ + fprintf(stderr,"** load_NBL_bricks: missing index list\n"); + return -1; + } + + prev = -1; /* use prev for previous sub-brick */ + for( c = 0; c < NBL->nbricks; c++ ){ + isrc = slist[c]; /* this is original brick index (c is new one) */ + idest = sindex[c]; /* this is the destination index for this data */ + + /* if this sub-brick is not the previous, we must read from disk */ + if( isrc != prev ){ + + /* if we are not looking at the correct sub-brick, scan forward */ + if( fposn != (oposn + isrc*NBL->bsize) ){ + fposn = oposn + isrc*NBL->bsize; + if( znzseek(fp, fposn, SEEK_SET) < 0 ){ + fprintf(stderr,"** failed to locate brick %d in file '%s'\n", + isrc, nim->iname ? nim->iname : nim->fname); + return -1; + } + } + + /* only 10,000 lines later and we're actually reading something! */ + rv = nifti_read_buffer(fp, NBL->bricks[idest], NBL->bsize, nim); + if( rv != NBL->bsize ){ + fprintf(stderr,"** failed to read brick %d from file '%s'\n", + isrc, nim->iname ? nim->iname : nim->fname); + return -1; + } + fposn += NBL->bsize; + } else { + /* we have already read this sub-brick, just copy the previous one */ + /* note that this works because they are sorted */ + memcpy(NBL->bricks[idest], NBL->bricks[sindex[c-1]], NBL->bsize); + } + + prev = isrc; /* in any case, note the now previous sub-brick */ + } + + return 0; +} + + +/*---------------------------------------------------------------------- + * nifti_alloc_NBL_mem - allocate memory for bricks + * + * return 0 on success, -1 on failure + *----------------------------------------------------------------------*/ +static int nifti_alloc_NBL_mem(nifti_image * nim, int nbricks, + nifti_brick_list * nbl) +{ + int c; + + /* if nbricks is not specified, use the default */ + if( nbricks > 0 ) nbl->nbricks = nbricks; + else { /* I missed this one with the 1.17 change 02 Mar 2006 [rickr] */ + nbl->nbricks = 1; + for( c = 4; c <= nim->ndim; c++ ) + nbl->nbricks *= nim->dim[c]; + } + + nbl->bsize = nim->nx * nim->ny * nim->nz * nim->nbyper; /* bytes */ + nbl->bricks = (void **)malloc(nbl->nbricks * sizeof(void *)); + + if( ! nbl->bricks ){ + fprintf(stderr,"** NANM: failed to alloc %d void ptrs\n",nbricks); + return -1; + } + + for( c = 0; c < nbl->nbricks; c++ ){ + nbl->bricks[c] = (void *)malloc(nbl->bsize); + if( ! nbl->bricks[c] ){ + fprintf(stderr,"** NANM: failed to alloc %d bytes for brick %d\n", + nbl->bsize, c); + /* so free and clear everything before returning */ + while( c > 0 ){ + c--; + free(nbl->bricks[c]); + } + free(nbl->bricks); + nbl->bricks = NULL; + nbl->nbricks = nbl->bsize = 0; + return -1; + } + } + + if( g_opts.debug > 2 ) + fprintf(stderr,"+d NANM: alloc'd %d bricks of %d bytes for NBL\n", + nbl->nbricks, nbl->bsize); + + return 0; +} + + +/*---------------------------------------------------------------------- + * nifti_copynsort - copy int list, and sort with indices + * + * 1. duplicate the incoming list + * 2. create an sindex list, and init with 0..nbricks-1 + * 3. do a slow insertion sort on the small slist, along with sindex list + * 4. check results, just to be positive + * + * So slist is sorted, and sindex hold original positions. + * + * return 0 on success, -1 on failure + *----------------------------------------------------------------------*/ +static int nifti_copynsort(int nbricks, const int * blist, int ** slist, + int ** sindex) +{ + int * stmp, * itmp; /* for ease of typing/reading */ + int c1, c2, spos, tmp; + + *slist = (int *)malloc(nbricks * sizeof(int)); + *sindex = (int *)malloc(nbricks * sizeof(int)); + + if( !*slist || !*sindex ){ + fprintf(stderr,"** NCS: failed to alloc %d ints for sorting\n",nbricks); + if(*slist) free(*slist); /* maybe one succeeded */ + if(*sindex) free(*sindex); + return -1; + } + + /* init the lists */ + memcpy(*slist, blist, nbricks*sizeof(int)); + for( c1 = 0; c1 < nbricks; c1++ ) (*sindex)[c1] = c1; + + /* now actually sort slist */ + stmp = *slist; + itmp = *sindex; + for( c1 = 0; c1 < nbricks-1; c1++ ) { + /* find smallest value, init to current */ + spos = c1; + for( c2 = c1+1; c2 < nbricks; c2++ ) + if( stmp[c2] < stmp[spos] ) spos = c2; + if( spos != c1 ) /* swap: fine, don't maintain sub-order, see if I care */ + { + tmp = stmp[c1]; /* first swap the sorting values */ + stmp[c1] = stmp[spos]; + stmp[spos] = tmp; + + tmp = itmp[c1]; /* then swap the index values */ + itmp[c1] = itmp[spos]; + itmp[spos] = tmp; + } + } + + if( g_opts.debug > 2 ){ + fprintf(stderr, "+d sorted indexing list:\n"); + fprintf(stderr, " orig : "); + for( c1 = 0; c1 < nbricks; c1++ ) fprintf(stderr," %d",blist[c1]); + fprintf(stderr,"\n new : "); + for( c1 = 0; c1 < nbricks; c1++ ) fprintf(stderr," %d",stmp[c1]); + fprintf(stderr,"\n indices: "); + for( c1 = 0; c1 < nbricks; c1++ ) fprintf(stderr," %d",itmp[c1]); + fputc('\n', stderr); + } + + /* check the sort (why not? I've got time...) */ + for( c1 = 0; c1 < nbricks-1; c1++ ){ + if( (stmp[c1] > stmp[c1+1]) || (blist[itmp[c1]] != stmp[c1]) ){ + fprintf(stderr,"** sorting screw-up, way to go, rick!\n"); + free(stmp); free(itmp); *slist = NULL; *sindex = NULL; + return -1; + } + } + + if( g_opts.debug > 2 ) fprintf(stderr,"-d sorting is okay\n"); + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! valid_nifti_brick_list - check sub-brick list for image + * + * This function verifies that nbricks and blist are appropriate + * for use with this nim, based on the dimensions. + * + * \param nim nifti_image to check against + * \param nbricks number of brick indices in blist + * \param blist list of brick indices to check in nim + * \param disp_error if this flag is set, report errors to user + * + * \return 1 if valid, 0 if not +*//*--------------------------------------------------------------------*/ +int valid_nifti_brick_list(nifti_image * nim , int nbricks, + const int * blist, int disp_error) +{ + int c, nsubs; + + if( !nim ){ + if( disp_error || g_opts.debug > 0 ) + fprintf(stderr,"** valid_nifti_brick_list: missing nifti image\n"); + return 0; + } + + if( nbricks <= 0 || !blist ){ + if( disp_error || g_opts.debug > 1 ) + fprintf(stderr,"** valid_nifti_brick_list: no brick list to check\n"); + return 0; + } + + /* nsubs sub-brick is nt*nu*nv*nw */ + for( c = 4, nsubs = 1; c <= nim->dim[0]; c++ ) + nsubs *= nim->dim[c]; + + if( nsubs <= 0 ){ + fprintf(stderr,"** VNBL warning: bad dim list (%d,%d,%d,%d)\n", + nim->dim[4], nim->dim[5], nim->dim[6], nim->dim[7]); + return 0; + } + + for( c = 0; c < nbricks; c++ ) + if( (blist[c] < 0) || (blist[c] >= nsubs) ){ + if( disp_error || g_opts.debug > 1 ) + fprintf(stderr, + "-d ** bad sub-brick chooser %d (#%d), valid range is [0,%d]\n", + blist[c], c, nsubs-1); + return 0; + } + + return 1; /* all is well */ +} + +#if 0 +/* set any non-positive values to 1 */ +static int int_force_positive( int * list, int nel ) +{ + int c; + if( !list || nel < 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** int_force_positive: bad params (%p,%d)\n", + (void *)list,nel); + return -1; + } + for( c = 0; c < nel; c++ ) + if( list[c] <= 0 ) list[c] = 1; + return 0; +} +#endif +/* end of new nifti_image_read_bricks() functionality */ + +/*----------------------------------------------------------------------*/ +/*! display the orientation from the quaternian fields + * + * \param mesg if non-NULL, display this message first + * \param mat the matrix to convert to "nearest" orientation + * + * \return -1 if results cannot be determined, 0 if okay +*//*--------------------------------------------------------------------*/ +int nifti_disp_matrix_orient( const char * mesg, mat44 mat ) +{ + int i, j, k; + + if ( mesg ) fputs( mesg, stderr ); /* use stdout? */ + + nifti_mat44_to_orientation( mat, &i,&j,&k ); + if ( i <= 0 || j <= 0 || k <= 0 ) return -1; + + /* so we have good codes */ + fprintf(stderr, " i orientation = '%s'\n" + " j orientation = '%s'\n" + " k orientation = '%s'\n", + nifti_orientation_string(i), + nifti_orientation_string(j), + nifti_orientation_string(k) ); + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! duplicate the given string (alloc length+1) + * + * \return allocated pointer (or NULL on failure) +*//*--------------------------------------------------------------------*/ +char *nifti_strdup(const char *str) +{ + char *dup= (char *)malloc( strlen(str)+1 ); + if (dup) strcpy(dup,str); + return dup; +} + + +/*---------------------------------------------------------------------------*/ +/*! Return a pointer to a string holding the name of a NIFTI datatype. + + \param dt NIfTI-1 datatype + + \return pointer to static string holding the datatype name + + \warning Do not free() or modify this string! + It points to static storage. + + \sa NIFTI1_DATATYPES group in nifti1.h +*//*-------------------------------------------------------------------------*/ char *nifti_datatype_string( int dt ) { switch( dt ){ @@ -43,10 +1036,13 @@ return "**ILLEGAL**" ; } -/*---------------------------------------------------------------------------*/ -/* Determine if the datatype code dt is an integer type (1=YES, 0=NO). ------------------------------------------------------------------------------*/ - +/*----------------------------------------------------------------------*/ +/*! Determine if the datatype code dt is an integer type (1=YES, 0=NO). + + \return whether the given NIfTI-1 datatype code is valid + + \sa NIFTI1_DATATYPES group in nifti1.h +*//*--------------------------------------------------------------------*/ int nifti_is_inttype( int dt ) { switch( dt ){ @@ -72,10 +1068,17 @@ } /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI units type. - Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/*! Return a pointer to a string holding the name of a NIFTI units type. + + \param uu NIfTI-1 unit code + + \return pointer to static string for the given unit type + + \warning Do not free() or modify this string! + It points to static storage. + + \sa NIFTI1_UNITS group in nifti1.h +*//*-------------------------------------------------------------------------*/ char *nifti_units_string( int uu ) { switch( uu ){ @@ -87,15 +1090,23 @@ case NIFTI_UNITS_USEC: return "us" ; case NIFTI_UNITS_HZ: return "Hz" ; case NIFTI_UNITS_PPM: return "ppm" ; + case NIFTI_UNITS_RADS: return "rad/s" ; } return "Unknown" ; } /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI transform type. - Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/*! Return a pointer to a string holding the name of a NIFTI transform type. + + \param xx NIfTI-1 xform code + + \return pointer to static string describing xform code + + \warning Do not free() or modify this string! + It points to static storage. + + \sa NIFTI1_XFORM_CODES group in nifti1.h +*//*-------------------------------------------------------------------------*/ char *nifti_xform_string( int xx ) { switch( xx ){ @@ -108,10 +1119,17 @@ } /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI intent type. - Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/*! Return a pointer to a string holding the name of a NIFTI intent type. + + \param ii NIfTI-1 intent code + + \return pointer to static string describing code + + \warning Do not free() or modify this string! + It points to static storage. + + \sa NIFTI1_INTENT_CODES group in nifti1.h +*//*-------------------------------------------------------------------------*/ char *nifti_intent_string( int ii ) { switch( ii ){ @@ -137,6 +1155,9 @@ case NIFTI_INTENT_EXTVAL: return "Extreme Value distribution" ; case NIFTI_INTENT_PVAL: return "P-value" ; + case NIFTI_INTENT_LOGPVAL: return "Log P-value" ; + case NIFTI_INTENT_LOG10PVAL: return "Log10 P-value" ; + case NIFTI_INTENT_ESTIMATE: return "Estimate" ; case NIFTI_INTENT_LABEL: return "Label index" ; case NIFTI_INTENT_NEURONAME: return "NeuroNames index" ; @@ -147,31 +1168,49 @@ case NIFTI_INTENT_POINTSET: return "Pointset" ; case NIFTI_INTENT_TRIANGLE: return "Triangle" ; case NIFTI_INTENT_QUATERNION: return "Quaternion" ; + + case NIFTI_INTENT_DIMLESS: return "Dimensionless number" ; } return "Unknown" ; } /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI slice_code - string. Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/*! Return a pointer to a string holding the name of a NIFTI slice_code. + + \param ss NIfTI-1 slice order code + + \return pointer to static string describing code + + \warning Do not free() or modify this string! + It points to static storage. + + \sa NIFTI1_SLICE_ORDER group in nifti1.h +*//*-------------------------------------------------------------------------*/ char *nifti_slice_string( int ss ) { switch( ss ){ - case NIFTI_SLICE_SEQ_INC: return "sequential_increasing" ; - case NIFTI_SLICE_SEQ_DEC: return "sequential_decreasing" ; - case NIFTI_SLICE_ALT_INC: return "alternating_increasing" ; - case NIFTI_SLICE_ALT_DEC: return "alternating_decreasing" ; + case NIFTI_SLICE_SEQ_INC: return "sequential_increasing" ; + case NIFTI_SLICE_SEQ_DEC: return "sequential_decreasing" ; + case NIFTI_SLICE_ALT_INC: return "alternating_increasing" ; + case NIFTI_SLICE_ALT_DEC: return "alternating_decreasing" ; + case NIFTI_SLICE_ALT_INC2: return "alternating_increasing_2" ; + case NIFTI_SLICE_ALT_DEC2: return "alternating_decreasing_2" ; } return "Unknown" ; } /*---------------------------------------------------------------------------*/ -/* Return a pointer to a string holding the name of a NIFTI orientation - string. Don't free() or modify this string! It points to static storage. ------------------------------------------------------------------------------*/ - +/*! Return a pointer to a string holding the name of a NIFTI orientation. + + \param ii orientation code + + \return pointer to static string holding the orientation information + + \warning Do not free() or modify the return string! + It points to static storage. + + \sa NIFTI_L2R in nifti1_io.h +*//*-------------------------------------------------------------------------*/ char *nifti_orientation_string( int ii ) { switch( ii ){ @@ -186,10 +1225,18 @@ } /*--------------------------------------------------------------------------*/ -/* Given a datatype code, set number of bytes per voxel and the swapsize. - The swapsize is set to 0 if this datatype doesn't ever need swapping. -----------------------------------------------------------------------------*/ - +/*! Given a datatype code, set number of bytes per voxel and the swapsize. + + \param datatype nifti1 datatype code + \param nbyper pointer to return value: number of bytes per voxel + \param swapsize pointer to return value: size of swap blocks + + \return appropriate values at nbyper and swapsize + + The swapsize is set to 0 if this datatype doesn't ever need swapping. + + \sa NIFTI1_DATATYPES in nifti1.h +*//*------------------------------------------------------------------------*/ void nifti_datatype_sizes( int datatype , int *nbyper, int *swapsize ) { int nb=0, ss=0 ; @@ -223,20 +1270,28 @@ } /*---------------------------------------------------------------------------*/ -/* Given the quaternion parameters (etc.), compute a transformation matrix. +/*! Given the quaternion parameters (etc.), compute a transformation matrix. + See comments in nifti1.h for details. - qb,qc,qd = quaternion parameters - qx,qy,qz = offset parameters - dx,dy,dz = grid stepsizes (non-negative inputs are set to 1.0) - qfac = sign of dz step (< 0 is negative; >= 0 is positive) + + <pre> If qx=qy=qz=0, dx=dy=dz=1, then the output is a rotation matrix. For qfac >= 0, the rotation is proper. For qfac < 0, the rotation is improper. ------------------------------------------------------------------------------*/ - -mat44 quatern_to_mat44( float qb, float qc, float qd, - float qx, float qy, float qz, - float dx, float dy, float dz, float qfac ) + </pre> + + \see "QUATERNION REPRESENTATION OF ROTATION MATRIX" in nifti1.h + \see nifti_mat44_to_quatern, nifti_make_orthog_mat44, + nifti_mat44_to_orientation + +*//*-------------------------------------------------------------------------*/ +mat44 nifti_quatern_to_mat44( float qb, float qc, float qd, + float qx, float qy, float qz, + float dx, float dy, float dz, float qfac ) { mat44 R ; double a,b=qb,c=qc,d=qd , xd,yd,zd ; @@ -282,8 +1337,9 @@ } /*---------------------------------------------------------------------------*/ -/* Given the 3x4 upper corner of the matrix R, compute the quaternion - parameters that fit it. See comments in nifti1.h for details. +/*! Given the 3x4 upper corner of the matrix R, compute the quaternion + parameters that fit it. + - Any NULL pointer on input won't get assigned (e.g., if you don't want dx,dy,dz, just pass NULL in for those pointers). - If the 3 input matrix columns are NOT orthogonal, they will be @@ -291,19 +1347,23 @@ the polar decomposition to find the orthogonal matrix closest to the column-normalized input matrix. - However, if the 3 input matrix columns are NOT orthogonal, then - the matrix produced by quatern_to_mat44 WILL have orthogonal + the matrix produced by nifti_quatern_to_mat44 WILL have orthogonal columns, so it won't be the same as the matrix input here. This "feature" is because the NIFTI 'qform' transform is deliberately not fully general -- it is intended to model a volume with perpendicular axes. - If the 3 input matrix columns are not even linearly independent, you'll just have to take your luck, won't you? ------------------------------------------------------------------------------*/ - -void mat44_to_quatern( mat44 R , - float *qb, float *qc, float *qd, - float *qx, float *qy, float *qz, - float *dx, float *dy, float *dz, float *qfac ) + + \see "QUATERNION REPRESENTATION OF ROTATION MATRIX" in nifti1.h + + \see nifti_quatern_to_mat44, nifti_make_orthog_mat44, + nifti_mat44_to_orientation +*//*-------------------------------------------------------------------------*/ +void nifti_mat44_to_quatern( mat44 R , + float *qb, float *qc, float *qd, + float *qx, float *qy, float *qz, + float *dx, float *dy, float *dz, float *qfac ) { double r11,r12,r13 , r21,r22,r23 , r31,r32,r33 ; double xd,yd,zd , a,b,c,d ; @@ -357,7 +1417,7 @@ Q.m[1][0] = r21 ; Q.m[1][1] = r22 ; Q.m[1][2] = r23 ; Q.m[2][0] = r31 ; Q.m[2][1] = r32 ; Q.m[2][2] = r33 ; - P = mat33_polar(Q) ; /* P is orthog matrix closest to Q */ + P = nifti_mat33_polar(Q) ; /* P is orthog matrix closest to Q */ r11 = P.m[0][0] ; r12 = P.m[0][1] ; r13 = P.m[0][2] ; /* unload */ r21 = P.m[1][0] ; r22 = P.m[1][1] ; r23 = P.m[1][2] ; @@ -416,18 +1476,26 @@ } /*---------------------------------------------------------------------------*/ -/* Compute the inverse of a bordered 4x4 matrix. - Some numerical code fragments were generated by Maple 8. - If a singular matrix is input, the output matrix will be all zero. - You can check for this by examining the [3][3] element, which will - be 1.0 for the normal case and 0.0 for the bad case. ------------------------------------------------------------------------------*/ - -mat44 mat44_inverse( mat44 R ) +/*! Compute the inverse of a bordered 4x4 matrix. + + <pre> + - Some numerical code fragments were generated by Maple 8. + - If a singular matrix is input, the output matrix will be all zero. + - You can check for this by examining the [3][3] element, which will + be 1.0 for the normal case and 0.0 for the bad case. + + The input matrix should have the form: + [ r11 r12 r13 v1 ] + [ r21 r22 r23 v2 ] + [ r31 r32 r33 v3 ] + [ 0 0 0 1 ] + </pre> +*//*-------------------------------------------------------------------------*/ +mat44 nifti_mat44_inverse( mat44 R ) { double r11,r12,r13,r21,r22,r23,r31,r32,r33,v1,v2,v3 , deti ; mat44 Q ; - /** INPUT MATRIX IS: **/ + /* INPUT MATRIX IS: */ r11 = R.m[0][0]; r12 = R.m[0][1]; r13 = R.m[0][2]; /* [ r11 r12 r13 v1 ] */ r21 = R.m[1][0]; r22 = R.m[1][1]; r23 = R.m[1][2]; /* [ r21 r22 r23 v2 ] */ r31 = R.m[2][0]; r32 = R.m[2][1]; r33 = R.m[2][2]; /* [ r31 r32 r33 v3 ] */ @@ -463,10 +1531,11 @@ } /*---------------------------------------------------------------------------*/ -/* Input 9 floats and make an orthgonal mat44 out of them. - Each row is normalized, then mat33_polar() is used to orthogonalize them. - If row #3 (r31,r32,r33) is input as zero, then it will be taken to be - the cross product of rows #1 and #2. +/*! Input 9 floats and make an orthgonal mat44 out of them. + + Each row is normalized, then nifti_mat33_polar() is used to orthogonalize + them. If row #3 (r31,r32,r33) is input as zero, then it will be taken to + be the cross product of rows #1 and #2. This function can be used to create a rotation matrix for transforming an oblique volume to anatomical coordinates. For this application: @@ -474,6 +1543,7 @@ - row #2 (r21,r22,r23) is the direction vector along the image j-axis - row #3 (r31,r32,r33) is the direction vector along the slice direction (if available; otherwise enter it as 0's) + The first 2 rows can be taken from the DICOM attribute (0020,0037) "Image Orientation (Patient)". @@ -483,13 +1553,16 @@ - column #1 (R.m[0][0],R.m[1][0],R.m[2][0]) by delta-x - column #2 (R.m[0][1],R.m[1][1],R.m[2][1]) by delta-y - column #3 (R.m[0][2],R.m[1][2],R.m[2][2]) by delta-z - And by placing the center (x,y,z) coordinates of voxel (0,0,0) into + + and by then placing the center (x,y,z) coordinates of voxel (0,0,0) into the column #4 (R.m[0][3],R.m[1][3],R.m[2][3]). ------------------------------------------------------------------------------*/ - -mat44 make_orthog_mat44( float r11, float r12, float r13 , - float r21, float r22, float r23 , - float r31, float r32, float r33 ) + + \sa nifti_quatern_to_mat44, nifti_mat44_to_quatern, + nifti_mat44_to_orientation +*//*-------------------------------------------------------------------------*/ +mat44 nifti_make_orthog_mat44( float r11, float r12, float r13 , + float r21, float r22, float r23 , + float r31, float r32, float r33 ) { mat44 R ; mat33 Q , P ; @@ -533,7 +1606,7 @@ Q.m[2][2] = Q.m[0][0]*Q.m[1][1] - Q.m[0][1]*Q.m[1][0] ; } - P = mat33_polar(Q) ; /* P is orthog matrix closest to Q */ + P = nifti_mat33_polar(Q) ; /* P is orthog matrix closest to Q */ R.m[0][0] = P.m[0][0] ; R.m[0][1] = P.m[0][1] ; R.m[0][2] = P.m[0][2] ; R.m[1][0] = P.m[1][0] ; R.m[1][1] = P.m[1][1] ; R.m[1][2] = P.m[1][2] ; @@ -542,13 +1615,14 @@ R.m[0][3] = R.m[1][3] = R.m[2][3] = 0.0 ; return R ; } -/*---------------------------------------------------------------------------*/ - -mat33 mat33_inverse( mat33 R ) /* inverse of 3x3 matrix */ +/*----------------------------------------------------------------------*/ +/*! compute the inverse of a 3x3 matrix +*//*--------------------------------------------------------------------*/ +mat33 nifti_mat33_inverse( mat33 R ) /* inverse of 3x3 matrix */ { double r11,r12,r13,r21,r22,r23,r31,r32,r33 , deti ; mat33 Q ; - /** INPUT MATRIX: **/ + /* INPUT MATRIX: */ r11 = R.m[0][0]; r12 = R.m[0][1]; r13 = R.m[0][2]; /* [ r11 r12 r13 ] */ r21 = R.m[1][0]; r22 = R.m[1][1]; r23 = R.m[1][2]; /* [ r21 r22 r23 ] */ r31 = R.m[2][0]; r32 = R.m[2][1]; r33 = R.m[2][2]; /* [ r31 r32 r33 ] */ @@ -573,12 +1647,13 @@ return Q ; } -/*---------------------------------------------------------------------------*/ - -float mat33_determ( mat33 R ) /* determinant of 3x3 matrix */ +/*----------------------------------------------------------------------*/ +/*! compute the determinant of a 3x3 matrix +*//*--------------------------------------------------------------------*/ +float nifti_mat33_determ( mat33 R ) /* determinant of 3x3 matrix */ { double r11,r12,r13,r21,r22,r23,r31,r32,r33 ; - /** INPUT MATRIX: **/ + /* INPUT MATRIX: */ r11 = R.m[0][0]; r12 = R.m[0][1]; r13 = R.m[0][2]; /* [ r11 r12 r13 ] */ r21 = R.m[1][0]; r22 = R.m[1][1]; r23 = R.m[1][2]; /* [ r21 r22 r23 ] */ r31 = R.m[2][0]; r32 = R.m[2][1]; r33 = R.m[2][2]; /* [ r31 r32 r33 ] */ @@ -587,9 +1662,10 @@ +r21*r32*r13+r31*r12*r23-r31*r22*r13 ; } -/*---------------------------------------------------------------------------*/ - -float mat33_rownorm( mat33 A ) /* max row norm of 3x3 matrix */ +/*----------------------------------------------------------------------*/ +/*! compute the max row norm of a 3x3 matrix +*//*--------------------------------------------------------------------*/ +float nifti_mat33_rownorm( mat33 A ) /* max row norm of 3x3 matrix */ { float r1,r2,r3 ; @@ -601,9 +1677,10 @@ return r1 ; } -/*---------------------------------------------------------------------------*/ - -float mat33_colnorm( mat33 A ) /* max column norm of 3x3 matrix */ +/*----------------------------------------------------------------------*/ +/*! compute the max column norm of a 3x3 matrix +*//*--------------------------------------------------------------------*/ +float nifti_mat33_colnorm( mat33 A ) /* max column norm of 3x3 matrix */ { float r1,r2,r3 ; @@ -615,11 +1692,12 @@ return r1 ; } -/*---------------------------------------------------------------------------*/ - -mat33 mat33_mul( mat33 A , mat33 B ) /* multiply 2 3x3 matrices */ +/*----------------------------------------------------------------------*/ +/*! multiply 2 3x3 matrices +*//*--------------------------------------------------------------------*/ +mat33 nifti_mat33_mul( mat33 A , mat33 B ) /* multiply 2 3x3 matrices */ { - mat33 C ; int i,j ; float sum ; + mat33 C ; int i,j ; for( i=0 ; i < 3 ; i++ ) for( j=0 ; j < 3 ; j++ ) C.m[i][j] = A.m[i][0] * B.m[0][j] @@ -629,12 +1707,14 @@ } /*---------------------------------------------------------------------------*/ -/* Polar decomposition of a 3x3 matrix: finds the closest orthogonal matrix - to input A (in both the Frobenius and L2 norms). Algorithm is that from - NJ Higham, SIAM J Sci Stat Comput, 7:1160-1174. ------------------------------------------------------------------------------*/ - -mat33 mat33_polar( mat33 A ) +/*! polar decomposition of a 3x3 matrix + + This finds the closest orthogonal matrix to input A + (in both the Frobenius and L2 norms). + + Algorithm is that from NJ Higham, SIAM J Sci Stat Comput, 7:1160-1174. +*//*-------------------------------------------------------------------------*/ +mat33 nifti_mat33_polar( mat33 A ) { mat33 X , Y , Z ; float alp,bet,gam,gmi , dif=1.0 ; @@ -644,18 +1724,18 @@ /* force matrix to be nonsingular */ - gam = mat33_determ(X) ; + gam = nifti_mat33_determ(X) ; while( gam == 0.0 ){ /* perturb matrix */ - gam = 0.00001 * ( 0.001 + mat33_rownorm(X) ) ; + gam = 0.00001 * ( 0.001 + nifti_mat33_rownorm(X) ) ; X.m[0][0] += gam ; X.m[1][1] += gam ; X.m[2][2] += gam ; - gam = mat33_determ(X) ; + gam = nifti_mat33_determ(X) ; } while(1){ - Y = mat33_inverse(X) ; + Y = nifti_mat33_inverse(X) ; if( dif > 0.3 ){ /* far from convergence */ - alp = sqrt( mat33_rownorm(X) * mat33_colnorm(X) ) ; - bet = sqrt( mat33_rownorm(Y) * mat33_colnorm(Y) ) ; + alp = sqrt( nifti_mat33_rownorm(X) * nifti_mat33_colnorm(X) ) ; + bet = sqrt( nifti_mat33_rownorm(Y) * nifti_mat33_colnorm(Y) ) ; gam = sqrt( bet / alp ) ; gmi = 1.0 / gam ; } else { @@ -686,7 +1766,10 @@ } /*---------------------------------------------------------------------------*/ -/* Input: 4x4 matrix that transforms (i,j,k) indexes to (x,y,z) coordinates, +/*! compute the (closest) orientation from a 4x4 ijk->xyz tranformation matrix + + <pre> + Input: 4x4 matrix that transforms (i,j,k) indexes to (x,y,z) coordinates, where +x=Right, +y=Anterior, +z=Superior. (Only the upper-left 3x3 corner of R is used herein.) Output: 3 orientation codes that correspond to the closest "standard" @@ -699,13 +1782,18 @@ *icod = NIFTI_R2L (i axis is mostly Right to Left) *jcod = NIFTI_P2A (j axis is mostly Posterior to Anterior) *kcod = NIFTI_I2S (k axis is mostly Inferior to Superior) ------------------------------------------------------------------------------*/ - -void mat44_to_orientation( mat44 R , int *icod, int *jcod, int *kcod ) + </pre> + + \see "QUATERNION REPRESENTATION OF ROTATION MATRIX" in nifti1.h + + \see nifti_quatern_to_mat44, nifti_mat44_to_quatern, + nifti_make_orthog_mat44 +*//*-------------------------------------------------------------------------*/ +void nifti_mat44_to_orientation( mat44 R , int *icod, int *jcod, int *kcod ) { float xi,xj,xk , yi,yj,yk , zi,zj,zk , val,detQ,detP ; mat33 P , Q , M ; - int i,j,k,p,q,r , ibest,jbest,kbest,pbest,qbest,rbest ; + int i,j,k=0,p,q,r , ibest,jbest,kbest,pbest,qbest,rbest ; float vbest ; if( icod == NULL || jcod == NULL || kcod == NULL ) return ; /* bad */ @@ -776,7 +1864,7 @@ /* at this point, Q is the rotation matrix from the (i,j,k) to (x,y,z) axes */ - detQ = mat33_determ( Q ) ; + detQ = nifti_mat33_determ( Q ) ; if( detQ == 0.0 ) return ; /* shouldn't happen unless user is a DUFIS */ /* Build and test all possible +1/-1 coordinate permutation matrices P; @@ -799,9 +1887,9 @@ for( q=-1 ; q <= 1 ; q+=2 ){ /* and go into rows #1,2,3 */ for( r=-1 ; r <= 1 ; r+=2 ){ P.m[0][i-1] = p ; P.m[1][j-1] = q ; P.m[2][k-1] = r ; - detP = mat33_determ(P) ; /* sign of permutation */ + detP = nifti_mat33_determ(P) ; /* sign of permutation */ if( detP * detQ <= 0.0 ) continue ; /* doesn't match sign of Q */ - M = mat33_mul(P,Q) ; + M = nifti_mat33_mul(P,Q) ; /* angle of M rotation = 2.0*acos(0.5*sqrt(1.0+trace(M))) */ /* we want largest trace(M) == smallest angle == M nearest to I */ @@ -868,7 +1956,10 @@ typedef struct { unsigned char a,b ; } twobytes ; -void swap_2bytes( int n , void *ar ) /* 2 bytes at a time */ +/*----------------------------------------------------------------------*/ +/*! swap each byte pair from the given list of n pairs +*//*--------------------------------------------------------------------*/ +void nifti_swap_2bytes( int n , void *ar ) /* 2 bytes at a time */ { register int ii ; register twobytes *tb = (twobytes *)ar ; @@ -884,7 +1975,10 @@ typedef struct { unsigned char a,b,c,d ; } fourbytes ; -void swap_4bytes( int n , void *ar ) /* 4 bytes at a time */ +/*----------------------------------------------------------------------*/ +/*! swap 4 bytes at a time from the given list of n sets of 4 bytes +*//*--------------------------------------------------------------------*/ +void nifti_swap_4bytes( int n , void *ar ) /* 4 bytes at a time */ { register int ii ; register fourbytes *tb = (fourbytes *)ar ; @@ -901,7 +1995,10 @@ typedef struct { unsigned char a,b,c,d , D,C,B,A ; } eightbytes ; -void swap_8bytes( int n , void *ar ) /* 8 bytes at a time */ +/*----------------------------------------------------------------------*/ +/*! swap 8 bytes at a time from the given list of n sets of 8 bytes +*//*--------------------------------------------------------------------*/ +void nifti_swap_8bytes( int n , void *ar ) /* 8 bytes at a time */ { register int ii ; register eightbytes *tb = (eightbytes *)ar ; @@ -921,7 +2018,10 @@ typedef struct { unsigned char a,b,c,d,e,f,g,h , H,G,F,E,D,C,B,A ; } sixteenbytes ; -void swap_16bytes( int n , void *ar ) /* 16 bytes at a time */ +/*----------------------------------------------------------------------*/ +/*! swap 16 bytes at a time from the given list of n sets of 16 bytes +*//*--------------------------------------------------------------------*/ +void nifti_swap_16bytes( int n , void *ar ) /* 16 bytes at a time */ { register int ii ; register sixteenbytes *tb = (sixteenbytes *)ar ; @@ -943,24 +2043,28 @@ /*---------------------------------------------------------------------------*/ -void swap_Nbytes( int n , int siz , void *ar ) /* subsuming case */ +/*----------------------------------------------------------------------*/ +/*! based on siz, call the appropriate nifti_swap_Nbytes() function +*//*--------------------------------------------------------------------*/ +void nifti_swap_Nbytes( int n , int siz , void *ar ) /* subsuming case */ { switch( siz ){ - case 2: swap_2bytes ( n , ar ) ; break ; - case 4: swap_4bytes ( n , ar ) ; break ; - case 8: swap_8bytes ( n , ar ) ; break ; - case 16: swap_16bytes( n , ar ) ; break ; + case 2: nifti_swap_2bytes ( n , ar ) ; break ; + case 4: nifti_swap_4bytes ( n , ar ) ; break ; + case 8: nifti_swap_8bytes ( n , ar ) ; break ; + case 16: nifti_swap_16bytes( n , ar ) ; break ; } return ; } -/*---------------------------------------------------------------*/ -/* Byte swap NIFTI-1 file header in various places and ways. - If is_nifti is nonzero, will also swap the NIFTI-specific - components of the header; otherwise, only the components - common to NIFTI and ANALYZE will be swapped. ----------------------------------------------------------------- */ - + +/*-------------------------------------------------------------------------*/ +/*! Byte swap NIFTI-1 file header in various places and ways. + + If is_nifti is nonzero, will also swap the NIFTI-specific + components of the header; otherwise, only the components + common to NIFTI and ANALYZE will be swapped. +*//*---------------------------------------------------------------------- */ void swap_nifti_header( struct nifti_1_header *h , int is_nifti ) { @@ -974,8 +2078,9 @@ /* this stuff is always present, for ANALYZE and NIFTI */ - swap_2bytes( 8 , h->dim ) ; - swap_4bytes( 8 , h->pixdim ) ; + swap_4(h->sizeof_hdr) ; + nifti_swap_2bytes( 8 , h->dim ) ; + nifti_swap_4bytes( 8 , h->pixdim ) ; swap_2(h->datatype) ; swap_2(h->bitpix) ; @@ -985,21 +2090,24 @@ /* this stuff is NIFTI specific */ if( is_nifti ){ + swap_4(h->intent_p1); swap_4(h->intent_p2); swap_4(h->intent_p3); + swap_2(h->intent_code); + + swap_2(h->slice_start); swap_2(h->slice_end); + swap_4(h->scl_slope); swap_4(h->scl_inter); + swap_4(h->slice_duration); swap_4(h->toffset); + swap_2(h->qform_code); swap_2(h->sform_code); swap_4(h->quatern_b); swap_4(h->quatern_c); swap_4(h->quatern_d); swap_4(h->qoffset_x); swap_4(h->qoffset_y); swap_4(h->qoffset_z); - swap_4(h->intent_p1); swap_4(h->intent_p2); swap_4(h->intent_p3); - swap_4(h->scl_slope); swap_4(h->scl_inter); - swap_4bytes(4,h->srow_x); - swap_4bytes(4,h->srow_y); - swap_4bytes(4,h->srow_z); - swap_2(h->intent_code); swap_4(h->toffset); - swap_2(h->slice_start); swap_2(h->slice_end); - swap_4(h->slice_duration); + nifti_swap_4bytes(4,h->srow_x); + nifti_swap_4bytes(4,h->srow_y); + nifti_swap_4bytes(4,h->srow_z); } return ; } + #define USE_STAT #ifdef USE_STAT /*---------------------------------------------------------------------------*/ @@ -1009,57 +2117,761 @@ #include <sys/types.h> #include <sys/stat.h> -unsigned int get_filesize( char *pathname ) +/*---------------------------------------------------------------------------*/ +/*! return the size of a file, in bytes + + \return size of file on success, -1 on error or no file + + changed to return int, -1 means no file or error 20 Dec 2004 [rickr] +*//*-------------------------------------------------------------------------*/ +int nifti_get_filesize( const char *pathname ) { struct stat buf ; int ii ; - if( pathname == NULL || *pathname == '\0' ) return 0 ; - ii = stat( pathname , &buf ); if( ii != 0 ) return 0 ; + if( pathname == NULL || *pathname == '\0' ) return -1 ; + ii = stat( pathname , &buf ); if( ii != 0 ) return -1 ; return (unsigned int)buf.st_size ; } #else /*---------- non-Unix version of the above, less efficient -----------*/ -unsigned int get_filesize( char *pathname ) +int nifti_get_filesize( const char *pathname ) { - FILE *fp ; unsigned int len ; - - if( pathname == NULL || *pathname == '\0' ) return 0 ; - fp = fopen(pathname,"rb"); if( fp == NULL ) return 0 ; - fseek(fp,0L,SEEK_END) ; len = (unsigned int)ftell(fp) ; - fclose(fp) ; return len ; + znzFile fp ; int len ; + + if( pathname == NULL || *pathname == '\0' ) return -1 ; + fp = znzopen(pathname,"rb",0); if( znz_isnull(fp) ) return -1 ; + znzseek(fp,0L,SEEK_END) ; len = znztell(fp) ; + znzclose(fp) ; return len ; } #endif /* USE_STAT */ + +/*----------------------------------------------------------------------*/ +/*! return the total volume size, in bytes + + This is computed as nvox * nbyper. +*//*--------------------------------------------------------------------*/ +size_t nifti_get_volsize(const nifti_image *nim) +{ + return (size_t)(nim->nbyper) * (size_t)(nim->nvox) ; /* total bytes */ +} + + /*--------------------------------------------------------------------------*/ -/* Determine if this is a NIFTI-formatted file. - - returns 0 if file looks like ANALYZE 7.5 [checks sizeof_hdr field == 348] - - returns 1 if file marked as NIFTI (header+data in 1 file) - - returns 2 if file marked as NIFTI (header+data in 2 files) - - returns -1 if it can't tell, file doesn't exist, etc. -----------------------------------------------------------------------------*/ - -int is_nifti_file( char *hname ) +/* Support functions for filenames in read and write + - allows for gzipped files +*/ + + +/*----------------------------------------------------------------------*/ +/*! simple check for file existence + + \return 1 on existence, 0 otherwise +*//*--------------------------------------------------------------------*/ +int nifti_fileexists(const char* fname) +{ + znzFile fp; + fp = znzopen( fname , "rb" , 1 ) ; + if( !znz_isnull(fp) ) { znzclose(fp); return 1; } + return 0; /* fp is NULL */ +} + +/*----------------------------------------------------------------------*/ +/*! return whether the filename is valid + + The name is considered valid if the file basename has length greater than + zero, AND one of the valid nifti extensions is provided. + fname input | return | + =============================== + "myimage" | 0 | + "myimage.tif" | 0 | + "myimage.tif.gz" | 0 | + "myimage.nii" | 1 | + ".nii" | 0 | + ".myhiddenimage" | 0 | + ".myhiddenimage.nii" | 1 | +*//*--------------------------------------------------------------------*/ +int nifti_is_complete_filename(const char* fname) +{ + char * ext; + + /* check input file(s) for sanity */ + if( fname == NULL || *fname == '\0' ){ + if ( g_opts.debug > 1 ) + fprintf(stderr,"-- empty filename in nifti_validfilename()\n"); + return 0; + } + + ext = nifti_find_file_extension(fname); + if ( ext == NULL ) { /*Invalid extension given */ + if ( g_opts.debug > 0 ) + fprintf(stderr,"-- no nifti valid extension for filename '%s'\n", fname); + return 0; + } + + if ( ext && ext == fname ) { /* then no filename prefix */ + if ( g_opts.debug > 0 ) + fprintf(stderr,"-- no prefix for filename '%s'\n", fname); + return 0; + } + return 1; +} + +/*----------------------------------------------------------------------*/ +/*! return whether the filename is valid + + The name is considered valid if its length is positive, excluding + any nifti filename extension. + fname input | return | result of nifti_makebasename + ==================================================================== + "myimage" | 1 | "myimage" + "myimage.tif" | 1 | "myimage.tif" + "myimage.tif.gz" | 1 | "myimage.tif" + "myimage.nii" | 1 | "myimage" + ".nii" | 0 | <ERROR - basename has zero length> + ".myhiddenimage" | 1 | ".myhiddenimage" + ".myhiddenimage.nii | 1 | ".myhiddenimage" +*//*--------------------------------------------------------------------*/ +int nifti_validfilename(const char* fname) +{ + char * ext; + + /* check input file(s) for sanity */ + if( fname == NULL || *fname == '\0' ){ + if ( g_opts.debug > 1 ) + fprintf(stderr,"-- empty filename in nifti_validfilename()\n"); + return 0; + } + + ext = nifti_find_file_extension(fname); + + if ( ext && ext == fname ) { /* then no filename prefix */ + if ( g_opts.debug > 0 ) + fprintf(stderr,"-- no prefix for filename '%s'\n", fname); + return 0; + } + + return 1; +} + +/*----------------------------------------------------------------------*/ +/*! check the end of the filename for a valid nifti extension + + Valid extensions are currently .nii, .hdr, .img, .nia, + or any of them followed by .gz. Note that '.' is part of + the extension. + + \return a pointer to the extension (within the filename), or NULL +*//*--------------------------------------------------------------------*/ +char * nifti_find_file_extension( const char * name ) +{ + char * ext; + int len; + + if ( ! name ) return NULL; + + len = strlen(name); + if ( len < 4 ) return NULL; + + ext = (char *)name + len - 4; + + if ( (strcmp(ext, ".hdr") == 0) || (strcmp(ext, ".img") == 0) || + (strcmp(ext, ".nia") == 0) || (strcmp(ext, ".nii") == 0) ) + return ext; + +#ifdef HAVE_ZLIB + if ( len < 7 ) return NULL; + + ext = (char *)name + len - 7; + + if ( (strcmp(ext, ".hdr.gz") == 0) || (strcmp(ext, ".img.gz") == 0) || + (strcmp(ext, ".nia.gz") == 0) || (strcmp(ext, ".nii.gz") == 0) ) + return ext; +#endif + + if( g_opts.debug > 1 ) + fprintf(stderr,"** find_file_ext: failed for name '%s'\n", name); + + return NULL; +} + +/*----------------------------------------------------------------------*/ +/*! return whether the filename ends in ".gz" +*//*--------------------------------------------------------------------*/ +int nifti_is_gzfile(const char* fname) +{ + /* return true if the filename ends with .gz */ + if (fname == NULL) { return 0; } +#ifdef HAVE_ZLIB + { /* just so len doesn't generate compile warning */ + int len; + len = strlen(fname); + if (len < 3) return 0; /* so we don't search before the name */ + if (strcmp(fname + strlen(fname) - 3,".gz")==0) { return 1; } + } +#endif + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! duplicate the filename, while clearing any extension + + This allocates memory for basename which should eventually be freed. +*//*--------------------------------------------------------------------*/ +char * nifti_makebasename(const char* fname) +{ + char *basename, *ext; + + basename=nifti_strdup(fname); + + ext = nifti_find_file_extension(basename); + if ( ext ) *ext = '\0'; /* clear out extension */ + + return basename; /* in either case */ +} + +/*----------------------------------------------------------------------*/ +/*! set nifti's global debug level, for status reporting + + - 0 : quiet, nothing is printed to the terminal, but errors + - 1 : normal execution (the default) + - 2, 3 : more details +*//*--------------------------------------------------------------------*/ +void nifti_set_debug_level( int level ) +{ + g_opts.debug = level; +} + +/*----------------------------------------------------------------------*/ +/*! set nifti's global skip_blank_ext flag 5 Sep 2006 [rickr] + + explicitly set to 0 or 1 +*//*--------------------------------------------------------------------*/ +void nifti_set_skip_blank_ext( int skip ) +{ + g_opts.skip_blank_ext = skip ? 1 : 0; +} + +/*----------------------------------------------------------------------*/ +/*! check current directory for existing header file + + \return filename of header on success and NULL if no appropriate file + could be found + + NB: it allocates memory for hdrname which should be freed + when no longer required +*//*-------------------------------------------------------------------*/ +char * nifti_findhdrname(const char* fname) +{ + char *basename, *hdrname, *ext; + char elist[2][5] = { ".hdr", ".nii" }; + int efirst; + + /**- check input file(s) for sanity */ + if( !nifti_validfilename(fname) ) return NULL; + + basename = nifti_makebasename(fname); + if( !basename ) return NULL; /* only on string alloc failure */ + + /**- return filename if it has a valid extension and exists + (except if it is an .img file (and maybe .gz)) */ + ext = nifti_find_file_extension(fname); + if ( ext && nifti_fileexists(fname) ) { + if ( strncmp(ext,".img",4) != 0 ){ + hdrname = nifti_strdup(fname); + free(basename); + return hdrname; + } + } + + /* So the requested name is a basename, contains .img, or does not exist. */ + /* In any case, use basename. */ + + /**- if .img, look for .hdr, .hdr.gz, .nii, .nii.gz, in that order */ + /**- else, look for .nii, .nii.gz, .hdr, .hdr.gz, in that order */ + + /* if we get more extension choices, this could be a loop */ + + if ( ext && strncmp(ext,".img",4) == 0 ) efirst = 0; + else efirst = 1; + + hdrname = (char *)calloc(sizeof(char),strlen(basename)+8); + if( !hdrname ){ + fprintf(stderr,"** nifti_findhdrname: failed to alloc hdrname\n"); + free(basename); + return NULL; + } + + strcpy(hdrname,basename); + strcat(hdrname,elist[efirst]); + if (nifti_fileexists(hdrname)) { free(basename); return hdrname; } +#ifdef HAVE_ZLIB + strcat(hdrname,".gz"); + if (nifti_fileexists(hdrname)) { free(basename); return hdrname; } +#endif + + /* okay, try the other possibility */ + + efirst = 1 - efirst; + + strcpy(hdrname,basename); + strcat(hdrname,elist[efirst]); + if (nifti_fileexists(hdrname)) { free(basename); return hdrname; } +#ifdef HAVE_ZLIB + strcat(hdrname,".gz"); + if (nifti_fileexists(hdrname)) { free(basename); return hdrname; } +#endif + + /**- if nothing has been found, return NULL */ + free(basename); + free(hdrname); + return NULL; +} + + +/*------------------------------------------------------------------------*/ +/*! check current directory for existing image file + + \param fname filename to check for + \nifti_type nifti_type for dataset - this determines whether to + first check for ".nii" or ".img" (since both may exist) + + \return filename of data/img file on success and NULL if no appropriate + file could be found + + NB: it allocates memory for the image filename, which should be freed + when no longer required +*//*---------------------------------------------------------------------*/ +char * nifti_findimgname(const char* fname , int nifti_type) +{ + char *basename, *imgname, ext[2][5] = { ".nii", ".img" }; + int first; /* first extension to use */ + + /* check input file(s) for sanity */ + + if( !nifti_validfilename(fname) ) return NULL; + + basename = nifti_makebasename(fname); + imgname = (char *)calloc(sizeof(char),strlen(basename)+8); + if( !imgname ){ + fprintf(stderr,"** nifti_findimgname: failed to alloc imgname\n"); + free(basename); + return NULL; + } + + /**- test for .nii and .img (don't assume input type from image type) */ + /**- if nifti_type = 1, check for .nii first, else .img first */ + + /* if we get 3 or more extensions, can make a loop here... */ + + if (nifti_type == NIFTI_FTYPE_NIFTI1_1) first = 0; /* should match .nii */ + else first = 1; /* should match .img */ + + strcpy(imgname,basename); + strcat(imgname,ext[first]); + if (nifti_fileexists(imgname)) { free(basename); return imgname; } +#ifdef HAVE_ZLIB /* then also check for .gz */ + strcat(imgname,".gz"); + if (nifti_fileexists(imgname)) { free(basename); return imgname; } +#endif + + /* failed to find image file with expected extension, try the other */ + + strcpy(imgname,basename); + strcat(imgname,ext[1-first]); /* can do this with only 2 choices */ + if (nifti_fileexists(imgname)) { free(basename); return imgname; } +#ifdef HAVE_ZLIB /* then also check for .gz */ + strcat(imgname,".gz"); + if (nifti_fileexists(imgname)) { free(basename); return imgname; } +#endif + + /**- if nothing has been found, return NULL */ + free(basename); + free(imgname); + return NULL; +} + + +/*----------------------------------------------------------------------*/ +/*! creates a filename for storing the header, based on nifti_type + + \param prefix - this will be copied before the suffix is added + \param nifti_type - determines the extension, unless one is in prefix + \param check - check for existence (fail condition) + \param comp - add .gz for compressed name + + Note that if prefix provides a file suffix, nifti_type is not used. + + NB: this allocates memory which should be freed + + \sa nifti_set_filenames +*//*-------------------------------------------------------------------*/ +char * nifti_makehdrname(const char * prefix, int nifti_type, int check, + int comp) +{ + char * iname, * ext; + + if( !nifti_validfilename(prefix) ) return NULL; + + /* add space for extension, optional ".gz", and null char */ + iname = (char *)calloc(sizeof(char),strlen(prefix)+8); + if( !iname ){ fprintf(stderr,"** small malloc failure!\n"); return NULL; } + strcpy(iname, prefix); + + /* use any valid extension */ + if( (ext = nifti_find_file_extension(iname)) != NULL ){ + if( strncmp(ext,".img",4) == 0 ) + memcpy(ext,".hdr",4); /* then convert img name to hdr */ + } + /* otherwise, make one up an */ + else if( nifti_type == NIFTI_FTYPE_NIFTI1_1 ) strcat(iname, ".nii"); + else if( nifti_type == NIFTI_FTYPE_ASCII ) strcat(iname, ".nia"); + else strcat(iname, ".hdr"); + +#ifdef HAVE_ZLIB /* then also check for .gz */ + if( comp && (!ext || !strstr(iname,".gz")) ) strcat(iname,".gz"); +#endif + + /* check for existence failure */ + if( check && nifti_fileexists(iname) ){ + fprintf(stderr,"** failure: header file '%s' already exists\n",iname); + free(iname); + return NULL; + } + + if(g_opts.debug > 2) fprintf(stderr,"+d made header filename '%s'\n", iname); + + return iname; +} + + +/*----------------------------------------------------------------------*/ +/*! creates a filename for storing the image, based on nifti_type + + \param prefix - this will be copied before the suffix is added + \param nifti_type - determines the extension, unless provided by prefix + \param check - check for existence (fail condition) + \param comp - add .gz for compressed name + + Note that if prefix provides a file suffix, nifti_type is not used. + + NB: it allocates memory which should be freed + + \sa nifti_set_filenames +*//*-------------------------------------------------------------------*/ +char * nifti_makeimgname(const char * prefix, int nifti_type, int check, + int comp) +{ + char * iname, * ext; + + if( !nifti_validfilename(prefix) ) return NULL; + + /* add space for extension, optional ".gz", and null char */ + iname = (char *)calloc(sizeof(char),strlen(prefix)+8); + if( !iname ){ fprintf(stderr,"** small malloc failure!\n"); return NULL; } + strcpy(iname, prefix); + + /* use any valid extension */ + if( (ext = nifti_find_file_extension(iname)) != NULL ){ + if( strncmp(ext,".hdr",4) == 0 ) + memcpy(ext,".img",4); /* then convert hdr name to img */ + } + /* otherwise, make one up */ + else if( nifti_type == NIFTI_FTYPE_NIFTI1_1 ) strcat(iname, ".nii"); + else if( nifti_type == NIFTI_FTYPE_ASCII ) strcat(iname, ".nia"); + else strcat(iname, ".img"); + +#ifdef HAVE_ZLIB /* then also check for .gz */ + if( comp && (!ext || !strstr(iname,".gz")) ) strcat(iname,".gz"); +#endif + + /* check for existence failure */ + if( check && nifti_fileexists(iname) ){ + fprintf(stderr,"** failure: image file '%s' already exists\n",iname); + free(iname); + return NULL; + } + + if( g_opts.debug > 2 ) fprintf(stderr,"+d made image filename '%s'\n",iname); + + return iname; +} + + +/*----------------------------------------------------------------------*/ +/*! create and set new filenames, based on prefix and image type + + \param nim pointer to nifti_image in which to set filenames + \param prefix (required) prefix for output filenames + \param check check for previous existence of filename + (existence is an error condition) + \param set_byte_order flag to set nim->byteorder here + (this is probably a logical place to do so) + + \return 0 on successful update + + \warning this will free() any existing names and create new ones + + \sa nifti_makeimgname, nifti_makehdrname, nifti_type_and_names_match +*//*--------------------------------------------------------------------*/ +int nifti_set_filenames( nifti_image * nim, const char * prefix, int check, + int set_byte_order ) +{ + int comp = nifti_is_gzfile(prefix); + + if( !nim || !prefix ){ + fprintf(stderr,"** nifti_set_filenames, bad params %p, %p\n", + (void *)nim,prefix); + return -1; + } + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d modifying output filenames using prefix %s\n", prefix); + + if( nim->fname ) free(nim->fname); + if( nim->iname ) free(nim->iname); + nim->fname = nifti_makehdrname(prefix, nim->nifti_type, check, comp); + nim->iname = nifti_makeimgname(prefix, nim->nifti_type, check, comp); + if( !nim->fname || !nim->iname ){ + LNI_FERR("nifti_set_filename","failed to set prefix for",prefix); + return -1; + } + + if( set_byte_order ) nim->byteorder = nifti_short_order() ; + + if( nifti_set_type_from_names(nim) < 0 ) + return -1; + + if( g_opts.debug > 2 ) + fprintf(stderr,"+d have new filenames %s and %s\n",nim->fname,nim->iname); + + return 0; +} + + +/*--------------------------------------------------------------------------*/ +/*! check whether nifti_type matches fname and iname for the nifti_image + + - if type 0 or 2, expect .hdr/.img pair + - if type 1, expect .nii (and names must match) + + \param nim given nifti_image + \param show_warn if set, print a warning message for any mis-match + + \return + - 1 if the values seem to match + - 0 if there is a mis-match + - -1 if there is not sufficient information to create file(s) + + \sa NIFTI_FTYPE_* codes in nifti1_io.h + \sa nifti_set_type_from_names, is_valid_nifti_type +*//*------------------------------------------------------------------------*/ +int nifti_type_and_names_match( nifti_image * nim, int show_warn ) +{ + char func[] = "nifti_type_and_names_match"; + char * ext_h, * ext_i; /* header and image filename extensions */ + int errs = 0; /* error counter */ + + /* sanity checks */ + if( !nim ){ + if( show_warn ) fprintf(stderr,"** %s: missing nifti_image\n", func); + return -1; + } + if( !nim->fname ){ + if( show_warn ) fprintf(stderr,"** %s: missing header filename\n", func); + errs++; + } + if( !nim->iname ){ + if( show_warn ) fprintf(stderr,"** %s: missing image filename\n", func); + errs++; + } + if( !is_valid_nifti_type(nim->nifti_type) ){ + if( show_warn ) + fprintf(stderr,"** %s: bad nifti_type %d\n", func, nim->nifti_type); + errs++; + } + + if( errs ) return -1; /* then do not proceed */ + + /* get pointers to extensions */ + ext_h = nifti_find_file_extension( nim->fname ); + ext_i = nifti_find_file_extension( nim->iname ); + + /* check for filename extensions */ + if( !ext_h ){ + if( show_warn ) + fprintf(stderr,"-d missing NIFTI extension in header filename, %s\n", + nim->fname); + errs++; + } + if( !ext_i ){ + if( show_warn ) + fprintf(stderr,"-d missing NIFTI extension in image filename, %s\n", + nim->iname); + errs++; + } + + if( errs ) return 0; /* do not proceed, but this is just a mis-match */ + + /* general tests */ + if( nim->nifti_type == NIFTI_FTYPE_NIFTI1_1 ){ /* .nii */ + if( strncmp(ext_h,".nii",4) ) { + if( show_warn ) + fprintf(stderr, + "-d NIFTI_FTYPE 1, but no .nii extension in header filename, %s\n", + nim->fname); + errs++; + } + if( strncmp(ext_i,".nii",4) ) { + if( show_warn ) + fprintf(stderr, + "-d NIFTI_FTYPE 1, but no .nii extension in image filename, %s\n", + nim->iname); + errs++; + } + if( strcmp(nim->fname, nim->iname) != 0 ){ + if( show_warn ) + fprintf(stderr, + "-d NIFTI_FTYPE 1, but header and image filenames differ: %s, %s\n", + nim->fname, nim->iname); + errs++; + } + } + else if( (nim->nifti_type == NIFTI_FTYPE_NIFTI1_2) || /* .hdr/.img */ + (nim->nifti_type == NIFTI_FTYPE_ANALYZE) ) + { + if( strncmp(ext_h,".hdr",4) != 0 ){ + if( show_warn ) + fprintf(stderr,"-d no '.hdr' extension, but NIFTI type is %d, %s\n", + nim->nifti_type, nim->fname); + errs++; + } + if( strncmp(ext_i,".img",4) != 0 ){ + if( show_warn ) + fprintf(stderr,"-d no '.img' extension, but NIFTI type is %d, %s\n", + nim->nifti_type, nim->iname); + errs++; + } + } + /* ignore any other nifti_type */ + + return 1; +} + + +/*--------------------------------------------------------------------------*/ +/*! check whether the given type is on the "approved" list + + The code is valid if it is non-negative, and does not exceed + NIFTI_MAX_FTYPE. + + \return 1 if nifti_type is valid, 0 otherwise + \sa NIFTI_FTYPE_* codes in nifti1_io.h +*//*------------------------------------------------------------------------*/ +int is_valid_nifti_type( int nifti_type ) +{ + if( nifti_type >= NIFTI_FTYPE_ANALYZE && /* smallest type, 0 */ + nifti_type <= NIFTI_MAX_FTYPE ) + return 1; + return 0; +} + + +/*--------------------------------------------------------------------------*/ +/*! set the nifti_type field based on fname and iname + + Note that nifti_type is changed only when it does not match + the filenames. + + \return 0 on success, -1 on error + + \sa is_valid_nifti_type, nifti_type_and_names_match +*//*------------------------------------------------------------------------*/ +int nifti_set_type_from_names( nifti_image * nim ) +{ + /* error checking first */ + if( !nim ){ fprintf(stderr,"** NSTFN: no nifti_image\n"); return -1; } + + if( !nim->fname || !nim->iname ){ + fprintf(stderr,"** NSTFN: missing filename(s) fname @ %p, iname @ %p\n", + nim->fname, nim->iname); + return -1; + } + + if( ! nifti_validfilename ( nim->fname ) || + ! nifti_validfilename ( nim->iname ) || + ! nifti_find_file_extension( nim->fname ) || + ! nifti_find_file_extension( nim->iname ) + ) { + fprintf(stderr,"** NSTFN: invalid filename(s) fname='%s', iname='%s'\n", + nim->fname, nim->iname); + return -1; + } + + if( g_opts.debug > 2 ) + fprintf(stderr,"-d verify nifti_type from filenames: %d",nim->nifti_type); + + /* not too picky here, do what must be done, and then verify */ + if( strcmp(nim->fname, nim->iname) == 0 ) /* one file, type 1 */ + nim->nifti_type = NIFTI_FTYPE_NIFTI1_1; + else if( nim->nifti_type == NIFTI_FTYPE_NIFTI1_1 ) /* cannot be type 1 */ + nim->nifti_type = NIFTI_FTYPE_NIFTI1_2; + + if( g_opts.debug > 2 ) fprintf(stderr," -> %d\n",nim->nifti_type); + + if( g_opts.debug > 1 ) /* warn user about anything strange */ + nifti_type_and_names_match(nim, 1); + + if( is_valid_nifti_type(nim->nifti_type) ) return 0; /* success! */ + + fprintf(stderr,"** NSTFN: bad nifti_type %d, for '%s' and '%s'\n", + nim->nifti_type, nim->fname, nim->iname); + + return -1; +} + + +/*--------------------------------------------------------------------------*/ +/*! Determine if this is a NIFTI-formatted file. + + <pre> + \return 0 if file looks like ANALYZE 7.5 [checks sizeof_hdr field == 348] + 1 if file marked as NIFTI (header+data in 1 file) + 2 if file marked as NIFTI (header+data in 2 files) + -1 if it can't tell, file doesn't exist, etc. + </pre> +*//*------------------------------------------------------------------------*/ +int is_nifti_file( const char *hname ) { struct nifti_1_header nhdr ; - FILE *fp ; + znzFile fp ; int ii ; + char *tmpname; /* bad input name? */ - if( hname == NULL || *hname == '\0' ) return -1 ; + if( !nifti_validfilename(hname) ) return -1 ; /* open file */ - fp = fopen( hname , "rb" ) ; - if( fp == NULL ) return -1 ; /* bad open? */ + tmpname = nifti_findhdrname(hname); + if( tmpname == NULL ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** no header file found for '%s'\n",hname); + return -1; + } + fp = znzopen( tmpname , "rb" , nifti_is_gzfile(tmpname) ) ; + free(tmpname); + if (znz_isnull(fp)) return -1 ; /* bad open? */ /* read header, close file */ - ii = fread( &nhdr , 1 , sizeof(nhdr) , fp ) ; - fclose( fp ) ; - if( ii < sizeof(nhdr) ) return -1 ; /* bad read? */ + ii = znzread( &nhdr , 1 , sizeof(nhdr) , fp ) ; + znzclose( fp ) ; + if( ii < (int) sizeof(nhdr) ) return -1 ; /* bad read? */ /* check for NIFTI-ness */ @@ -1070,405 +2882,1514 @@ /* check for ANALYZE-ness (sizeof_hdr field == 348) */ ii = nhdr.sizeof_hdr ; - if( ii == sizeof(nhdr) ) return 0 ; /* matches */ + if( ii == (int)sizeof(nhdr) ) return 0 ; /* matches */ /* try byte-swapping header */ swap_4(ii) ; - if( ii == sizeof(nhdr) ) return 0 ; /* matches */ + if( ii == (int)sizeof(nhdr) ) return 0 ; /* matches */ return -1 ; /* not good */ } -/*--------------------------------------------------------------------------*/ -/* Read in a NIFTI-1 or ANALYZE-7.5 file (pair) into a nifti_image struct. - - Input is .hdr or .nii filename. - - Return value is NULL if something fails badly. - - If read_data parameter is nonzero, the image data will actually - be read in; otherwise, it will have to be read later - (e.g., using the nifti_image_load() function). - - The image data will be stored in whatever data format the - input data is; no scaling will be applied. - - DT_BINARY data is not supported! - - nifti_image_free() can be used to delete the returned struct, - when you are done with it. -----------------------------------------------------------------------------*/ +static int print_hex_vals( const char * data, int nbytes, FILE * fp ) +{ + int c; + + if ( !data || nbytes < 1 || !fp ) return -1; + + fputs("0x", fp); + for ( c = 0; c < nbytes; c++ ) + fprintf(fp, " %x", data[c]); + + return 0; +} + +/*----------------------------------------------------------------------*/ +/*! display the contents of the nifti_1_header (send to stdout) + + \param info if non-NULL, print this character string + \param hp pointer to nifti_1_header +*//*--------------------------------------------------------------------*/ +int disp_nifti_1_header( const char * info, const nifti_1_header * hp ) +{ + int c; + + fputs( "-------------------------------------------------------\n", stdout ); + if ( info ) fputs( info, stdout ); + if ( !hp ){ fputs(" ** no nifti_1_header to display!\n",stdout); return 1; } + + fprintf(stdout," nifti_1_header :\n" + " sizeof_hdr = %d\n" + " data_type[10] = ", hp->sizeof_hdr); + print_hex_vals(hp->data_type, 10, stdout); + fprintf(stdout, "\n" + " db_name[18] = "); + print_hex_vals(hp->db_name, 18, stdout); + fprintf(stdout, "\n" + " extents = %d\n" + " session_error = %d\n" + " regular = 0x%x\n" + " dim_info = 0x%x\n", + hp->extents, hp->session_error, hp->regular, hp->dim_info ); + fprintf(stdout, " dim[8] ="); + for ( c = 0; c < 8; c++ ) fprintf(stdout," %d", hp->dim[c]); + fprintf(stdout, "\n" + " intent_p1 = %f\n" + " intent_p2 = %f\n" + " intent_p3 = %f\n" + " intent_code = %d\n" + " datatype = %d\n" + " bitpix = %d\n" + " slice_start = %d\n" + " pixdim[8] =", + hp->intent_p1, hp->intent_p2, hp->intent_p3, hp->intent_code, + hp->datatype, hp->bitpix, hp->slice_start); + /* break pixdim over 2 lines */ + for ( c = 0; c < 4; c++ ) fprintf(stdout," %f", hp->pixdim[c]); + fprintf(stdout, "\n "); + for ( c = 4; c < 8; c++ ) fprintf(stdout," %f", hp->pixdim[c]); + fprintf(stdout, "\n" + " vox_offset = %f\n" + " scl_slope = %f\n" + " scl_inter = %f\n" + " slice_end = %d\n" + " slice_code = %d\n" + " xyzt_units = 0x%x\n" + " cal_max = %f\n" + " cal_min = %f\n" + " slice_duration = %f\n" + " toffset = %f\n" + " glmax = %d\n" + " glmin = %d\n", + hp->vox_offset, hp->scl_slope, hp->scl_inter, hp->slice_end, + hp->slice_code, hp->xyzt_units, hp->cal_max, hp->cal_min, + hp->slice_duration, hp->toffset, hp->glmax, hp->glmin); + fprintf(stdout, + " descrip = '%.80s'\n" + " aux_file = '%.24s'\n" + " qform_code = %d\n" + " sform_code = %d\n" + " quatern_b = %f\n" + " quatern_c = %f\n" + " quatern_d = %f\n" + " qoffset_x = %f\n" + " qoffset_y = %f\n" + " qoffset_z = %f\n" + " srow_x[4] = %f, %f, %f, %f\n" + " srow_y[4] = %f, %f, %f, %f\n" + " srow_z[4] = %f, %f, %f, %f\n" + " intent_name = '%-.16s'\n" + " magic = '%-.4s'\n", + hp->descrip, hp->aux_file, hp->qform_code, hp->sform_code, + hp->quatern_b, hp->quatern_c, hp->quatern_d, + hp->qoffset_x, hp->qoffset_y, hp->qoffset_z, + hp->srow_x[0], hp->srow_x[1], hp->srow_x[2], hp->srow_x[3], + hp->srow_y[0], hp->srow_y[1], hp->srow_y[2], hp->srow_y[3], + hp->srow_z[0], hp->srow_z[1], hp->srow_z[2], hp->srow_z[3], + hp->intent_name, hp->magic); + fputs( "-------------------------------------------------------\n", stdout ); + fflush(stdout); + + return 0; +} + #undef ERREX #define ERREX(msg) \ - do{ fprintf(stderr,"** ERROR: nifti_image_read(%s): %s\n", \ - (hname != NULL) ? hname : "(null)" , (msg) ) ; \ + do{ fprintf(stderr,"** ERROR: nifti_convert_nhdr2nim: %s\n", (msg) ) ; \ return NULL ; } while(0) -nifti_image *nifti_image_read( char *hname , int read_data ) +/*----------------------------------------------------------------------*/ +/*! convert a nifti_1_header into a nift1_image + + \return an allocated nifti_image, or NULL on failure +*//*--------------------------------------------------------------------*/ +nifti_image* nifti_convert_nhdr2nim(struct nifti_1_header nhdr, + const char * fname) { - struct nifti_1_header nhdr ; - nifti_image *nim ; - FILE *fp ; - int ii , doswap , hlen, ilen, ioff , ndim,nvox ; + int ii , doswap , ioff ; int is_nifti , is_onefile ; - short ss ; - char *iname=NULL , buf[16] ; - - /** check input file(s) for sanity **/ - - if( hname == NULL || *hname == '\0' ) ERREX("bad filename") ; - - hlen = strlen(hname) ; if( hlen < 5 ) ERREX("too short filename") ; - - /** open input file **/ - - fp = fopen( hname , "rb" ) ; - if( fp == NULL ) ERREX("can't open header file") ; - - /** test if header file starts with ASCII string "<nifti_image"; - if so, read the dataset that special way and return now - (this is NOT part of the NIFTI-1 standard!) **/ - - ii = fread( buf , 1 , 12 , fp ) ; - if( ii < 12 ){ fclose( fp ) ; ERREX("bad header read") ; } - rewind(fp) ; - buf[12] = '\0' ; - if( strcmp(buf,"<nifti_image") == 0 ){ /* have an ASCII header! */ - int slen = get_filesize( hname ) ; - char *sbuf ; - if( slen > 65530 ) slen = 65530 ; - sbuf = (char *)calloc(sizeof(char),slen+1) ; - fread( sbuf , 1 , slen , fp ) ; fclose( fp ) ; - nim = nifti_image_from_ascii( sbuf ) ; free( sbuf ) ; - if( nim == NULL ) ERREX("bad ASCII header read") ; - nim->nifti_type = 3 ; - nim->iname_offset = -1 ; - if( read_data ) nifti_image_load( nim ) ; - else nim->data = NULL ; - return nim ; - } - - /***** Normal case: read binary header *****/ - - ii = fread( &nhdr , 1 , sizeof(nhdr) , fp ) ; /* read the thing */ - fclose( fp ) ; /* close the file */ - if( ii < sizeof(nhdr) ) ERREX("bad binary header read") ; - - /** check if have to swap bytes **/ - - doswap = 0 ; /* swap data flag */ - ss = nhdr.dim[0] ; - if( ss != 0 ){ /* check dim[0] for good value */ - if( ss < 0 || ss > 7 ){ - swap_2(ss) ; - if( ss < 0 || ss > 7 ) ERREX("bad dim[0]") ; - doswap = 1 ; - } - } else { /* dim[0] == 0 is illegal, but does occur */ - ii = nhdr.sizeof_hdr ; /* so check sizeof_hdr field instead */ - if( ii != sizeof(nhdr) ){ - swap_4(ii) ; - if( ii != sizeof(nhdr) ) ERREX("bad sizeof_hdr") ; - doswap = 1 ; - } - } - - /** determine if this is a NIFTI-1 compliant header **/ + nifti_image *nim; + + nim = (nifti_image *) calloc( 1 , sizeof(nifti_image) ) ; + if( !nim ) ERREX("failed to allocate nifti image"); + + /* be explicit with pointers */ + nim->fname = NULL; + nim->iname = NULL; + nim->data = NULL; + + /**- check if we must swap bytes */ + + doswap = need_nhdr_swap(nhdr.dim[0], nhdr.sizeof_hdr); /* swap data flag */ + + if( doswap < 0 ){ + if( doswap == -1 ) ERREX("bad dim[0]") ; + ERREX("bad sizeof_hdr") ; /* else */ + } + + /**- determine if this is a NIFTI-1 compliant header */ is_nifti = NIFTI_VERSION(nhdr) ; - if( doswap ) swap_nifti_header( &nhdr , is_nifti ) ; + if( doswap ) { + if ( g_opts.debug > 3 ) disp_nifti_1_header("-d ni1 pre-swap: ", &nhdr); + swap_nifti_header( &nhdr , is_nifti ) ; + } + + if ( g_opts.debug > 2 ) disp_nifti_1_header("-d nhdr2nim : ", &nhdr); if( nhdr.datatype == DT_BINARY || nhdr.datatype == DT_UNKNOWN ) ERREX("bad datatype") ; if( nhdr.dim[1] <= 0 ) ERREX("bad dim[1]") ; - for( ii=2 ; ii <= 7 ; ii++ ) - if( nhdr.dim[ii] <= 0 ) nhdr.dim[ii] = 1 ; /* fix bad dim[] values */ - - /** get number of dimensions (ignoring dim[0] now) **/ - + /* fix bad dim[] values in the defined dimension range */ + for( ii=2 ; ii <= nhdr.dim[0] ; ii++ ) + if( nhdr.dim[ii] <= 0 ) nhdr.dim[ii] = 1 ; + + /* fix any remaining bad dim[] values, so garbage does not propagate */ + /* (only values 0 or 1 seem rational, otherwise set to arbirary 1) */ + for( ii=nhdr.dim[0]+1 ; ii <= 7 ; ii++ ) + if( nhdr.dim[ii] != 1 && nhdr.dim[ii] != 0) nhdr.dim[ii] = 1 ; + +#if 0 /* rely on dim[0], do not attempt to modify it 16 Nov 2005 [rickr] */ + + /**- get number of dimensions (ignoring dim[0] now) */ for( ii=7 ; ii >= 2 ; ii-- ) /* loop backwards until we */ if( nhdr.dim[ii] > 1 ) break ; /* find a dim bigger than 1 */ ndim = ii ; - - /* set bad grid spacings to 1.0 */ - - for( ii=1 ; ii <= 7 ; ii++ ){ +#endif + + /**- set bad grid spacings to 1.0 */ + + for( ii=1 ; ii <= nhdr.dim[0] ; ii++ ){ if( nhdr.pixdim[ii] == 0.0 || !IS_GOOD_FLOAT(nhdr.pixdim[ii]) ) nhdr.pixdim[ii] = 1.0 ; } - /** will read image data from file 'iname' starting at offset 'ioff' **/ + is_onefile = is_nifti && NIFTI_ONEFILE(nhdr) ; + + if( is_nifti ) nim->nifti_type = (is_onefile) ? NIFTI_FTYPE_NIFTI1_1 + : NIFTI_FTYPE_NIFTI1_2 ; + else nim->nifti_type = NIFTI_FTYPE_ANALYZE ; + + ii = nifti_short_order() ; + if( doswap ) nim->byteorder = REVERSE_ORDER(ii) ; + else nim->byteorder = ii ; + + + /**- set dimensions of data array */ + + nim->ndim = nim->dim[0] = nhdr.dim[0]; + nim->nx = nim->dim[1] = nhdr.dim[1]; + nim->ny = nim->dim[2] = nhdr.dim[2]; + nim->nz = nim->dim[3] = nhdr.dim[3]; + nim->nt = nim->dim[4] = nhdr.dim[4]; + nim->nu = nim->dim[5] = nhdr.dim[5]; + nim->nv = nim->dim[6] = nhdr.dim[6]; + nim->nw = nim->dim[7] = nhdr.dim[7]; + + for( ii=1, nim->nvox=1; ii <= nhdr.dim[0]; ii++ ) + nim->nvox *= nhdr.dim[ii]; + + /**- set the type of data in voxels and how many bytes per voxel */ + + nim->datatype = nhdr.datatype ; + + nifti_datatype_sizes( nim->datatype , &(nim->nbyper) , &(nim->swapsize) ) ; + if( nim->nbyper == 0 ){ free(nim); ERREX("bad datatype"); } + + /**- set the grid spacings */ + + nim->dx = nim->pixdim[1] = nhdr.pixdim[1] ; + nim->dy = nim->pixdim[2] = nhdr.pixdim[2] ; + nim->dz = nim->pixdim[3] = nhdr.pixdim[3] ; + nim->dt = nim->pixdim[4] = nhdr.pixdim[4] ; + nim->du = nim->pixdim[5] = nhdr.pixdim[5] ; + nim->dv = nim->pixdim[6] = nhdr.pixdim[6] ; + nim->dw = nim->pixdim[7] = nhdr.pixdim[7] ; + + /**- compute qto_xyz transformation from pixel indexes (i,j,k) to (x,y,z) */ + + if( !is_nifti || nhdr.qform_code <= 0 ){ + /**- if not nifti or qform_code <= 0, use grid spacing for qto_xyz */ + + nim->qto_xyz.m[0][0] = nim->dx ; /* grid spacings */ + nim->qto_xyz.m[1][1] = nim->dy ; /* along diagonal */ + nim->qto_xyz.m[2][2] = nim->dz ; + + /* off diagonal is zero */ + + nim->qto_xyz.m[0][1]=nim->qto_xyz.m[0][2]=nim->qto_xyz.m[0][3] = 0.0; + nim->qto_xyz.m[1][0]=nim->qto_xyz.m[1][2]=nim->qto_xyz.m[1][3] = 0.0; + nim->qto_xyz.m[2][0]=nim->qto_xyz.m[2][1]=nim->qto_xyz.m[2][3] = 0.0; + + /* last row is always [ 0 0 0 1 ] */ + + nim->qto_xyz.m[3][0]=nim->qto_xyz.m[3][1]=nim->qto_xyz.m[3][2] = 0.0; + nim->qto_xyz.m[3][3]= 1.0 ; + + nim->qform_code = NIFTI_XFORM_UNKNOWN ; + + if( g_opts.debug > 1 ) fprintf(stderr,"-d no qform provided\n"); + } else { + /**- else NIFTI: use the quaternion-specified transformation */ + + nim->quatern_b = FIXED_FLOAT( nhdr.quatern_b ) ; + nim->quatern_c = FIXED_FLOAT( nhdr.quatern_c ) ; + nim->quatern_d = FIXED_FLOAT( nhdr.quatern_d ) ; + + nim->qoffset_x = FIXED_FLOAT(nhdr.qoffset_x) ; + nim->qoffset_y = FIXED_FLOAT(nhdr.qoffset_y) ; + nim->qoffset_z = FIXED_FLOAT(nhdr.qoffset_z) ; + + nim->qfac = (nhdr.pixdim[0] < 0.0) ? -1.0 : 1.0 ; /* left-handedness? */ + + nim->qto_xyz = nifti_quatern_to_mat44( + nim->quatern_b, nim->quatern_c, nim->quatern_d, + nim->qoffset_x, nim->qoffset_y, nim->qoffset_z, + nim->dx , nim->dy , nim->dz , + nim->qfac ) ; + + nim->qform_code = nhdr.qform_code ; + + if( g_opts.debug > 1 ) + nifti_disp_matrix_orient("-d qform orientations:\n", nim->qto_xyz); + } + + /**- load inverse transformation (x,y,z) -> (i,j,k) */ + + nim->qto_ijk = nifti_mat44_inverse( nim->qto_xyz ) ; + + /**- load sto_xyz affine transformation, if present */ + + if( !is_nifti || nhdr.sform_code <= 0 ){ + /**- if not nifti or sform_code <= 0, then no sto transformation */ + + nim->sform_code = NIFTI_XFORM_UNKNOWN ; + + if( g_opts.debug > 1 ) fprintf(stderr,"-d no sform provided\n"); + + } else { + /**- else set the sto transformation from srow_*[] */ + + nim->sto_xyz.m[0][0] = nhdr.srow_x[0] ; + nim->sto_xyz.m[0][1] = nhdr.srow_x[1] ; + nim->sto_xyz.m[0][2] = nhdr.srow_x[2] ; + nim->sto_xyz.m[0][3] = nhdr.srow_x[3] ; + + nim->sto_xyz.m[1][0] = nhdr.srow_y[0] ; + nim->sto_xyz.m[1][1] = nhdr.srow_y[1] ; + nim->sto_xyz.m[1][2] = nhdr.srow_y[2] ; + nim->sto_xyz.m[1][3] = nhdr.srow_y[3] ; + + nim->sto_xyz.m[2][0] = nhdr.srow_z[0] ; + nim->sto_xyz.m[2][1] = nhdr.srow_z[1] ; + nim->sto_xyz.m[2][2] = nhdr.srow_z[2] ; + nim->sto_xyz.m[2][3] = nhdr.srow_z[3] ; + + /* last row is always [ 0 0 0 1 ] */ + + nim->sto_xyz.m[3][0]=nim->sto_xyz.m[3][1]=nim->sto_xyz.m[3][2] = 0.0; + nim->sto_xyz.m[3][3]= 1.0 ; + + nim->sto_ijk = nifti_mat44_inverse( nim->sto_xyz ) ; + + nim->sform_code = nhdr.sform_code ; + + if( g_opts.debug > 1 ) + nifti_disp_matrix_orient("-d sform orientations:\n", nim->sto_xyz); + } + + /**- set miscellaneous NIFTI stuff */ + + if( is_nifti ){ + nim->scl_slope = FIXED_FLOAT( nhdr.scl_slope ) ; + nim->scl_inter = FIXED_FLOAT( nhdr.scl_inter ) ; + + nim->intent_code = nhdr.intent_code ; + + nim->intent_p1 = FIXED_FLOAT( nhdr.intent_p1 ) ; + nim->intent_p2 = FIXED_FLOAT( nhdr.intent_p2 ) ; + nim->intent_p3 = FIXED_FLOAT( nhdr.intent_p3 ) ; + + nim->toffset = FIXED_FLOAT( nhdr.toffset ) ; + + memcpy(nim->intent_name,nhdr.intent_name,15); nim->intent_name[15] = '\0'; + + nim->xyz_units = XYZT_TO_SPACE(nhdr.xyzt_units) ; + nim->time_units = XYZT_TO_TIME (nhdr.xyzt_units) ; + + nim->freq_dim = DIM_INFO_TO_FREQ_DIM ( nhdr.dim_info ) ; + nim->phase_dim = DIM_INFO_TO_PHASE_DIM( nhdr.dim_info ) ; + nim->slice_dim = DIM_INFO_TO_SLICE_DIM( nhdr.dim_info ) ; + + nim->slice_code = nhdr.slice_code ; + nim->slice_start = nhdr.slice_start ; + nim->slice_end = nhdr.slice_end ; + nim->slice_duration = FIXED_FLOAT(nhdr.slice_duration) ; + } + + /**- set Miscellaneous ANALYZE stuff */ + + nim->cal_min = FIXED_FLOAT(nhdr.cal_min) ; + nim->cal_max = FIXED_FLOAT(nhdr.cal_max) ; + + memcpy(nim->descrip ,nhdr.descrip ,79) ; nim->descrip [79] = '\0' ; + memcpy(nim->aux_file,nhdr.aux_file,23) ; nim->aux_file[23] = '\0' ; + + /**- set ioff from vox_offset (but at least sizeof(header)) */ is_onefile = is_nifti && NIFTI_ONEFILE(nhdr) ; if( is_onefile ){ ioff = (int)nhdr.vox_offset ; - if( ioff < sizeof(nhdr) ) ioff = sizeof(nhdr) ; + if( ioff < (int) sizeof(nhdr) ) ioff = (int) sizeof(nhdr) ; } else { - ioff = 0 ; - } - - iname = strdup(hname) ; - if( !is_onefile ) strcpy(iname+hlen-4,".img") ; /* create .img filename */ - - ilen = get_filesize(iname) ; /* find size of image data file */ - - if( ilen <= ioff ){ free(iname) ; ERREX("bad data file") ; } - - /*=== create output image struct and start to set it up ===*/ - - nim = (nifti_image *) calloc( 1 , sizeof(nifti_image) ) ; - - if( is_nifti ) nim->nifti_type = (is_onefile) ? 1 : 2 ; - else nim->nifti_type = 0 ; - - ii = short_order() ; - if( doswap ) nim->byteorder = REVERSE_ORDER(ii) ; - else nim->byteorder = ii ; - - /** dimensions of data array **/ - - nim->ndim = nim->dim[0] = ndim ; - nim->nx = nim->dim[1] = nhdr.dim[1]; nvox = nim->nx; - nim->ny = nim->dim[2] = nhdr.dim[2]; nvox *= nim->ny; - nim->nz = nim->dim[3] = nhdr.dim[3]; nvox *= nim->nz; - nim->nt = nim->dim[4] = nhdr.dim[4]; nvox *= nim->nt; - nim->nu = nim->dim[5] = nhdr.dim[5]; nvox *= nim->nu; - nim->nv = nim->dim[6] = nhdr.dim[6]; nvox *= nim->nv; - nim->nw = nim->dim[7] = nhdr.dim[7]; nvox *= nim->nw; nim->nvox = nvox; - - /** type of data in voxels and how many bytes per voxel */ - - nim->datatype = nhdr.datatype ; - - nifti_datatype_sizes( nim->datatype , &(nim->nbyper) , &(nim->swapsize) ) ; - if( nim->nbyper == 0 ){ free(nim); free(iname); ERREX("bad datatype"); } - - /** grid spacings **/ - - nim->dx = nim->pixdim[1] = nhdr.pixdim[1] ; - nim->dy = nim->pixdim[2] = nhdr.pixdim[2] ; - nim->dz = nim->pixdim[3] = nhdr.pixdim[3] ; - nim->dt = nim->pixdim[4] = nhdr.pixdim[4] ; - nim->du = nim->pixdim[5] = nhdr.pixdim[5] ; - nim->dv = nim->pixdim[6] = nhdr.pixdim[6] ; - nim->dw = nim->pixdim[7] = nhdr.pixdim[7] ; - - /** compute qto_xyz transformation from pixel indexes (i,j,k) to (x,y,z) **/ - - if( !is_nifti || nhdr.qform_code <= 0 ){ /** default transformation **/ - - nim->qto_xyz.m[0][0] = nim->dx ; /* grid spacings */ - nim->qto_xyz.m[1][1] = nim->dy ; /* along diagonal */ - nim->qto_xyz.m[2][2] = nim->dz ; - - /* off diagonal is zero */ - - nim->qto_xyz.m[0][1]=nim->qto_xyz.m[0][2]=nim->qto_xyz.m[0][3] = 0.0; - nim->qto_xyz.m[1][0]=nim->qto_xyz.m[1][2]=nim->qto_xyz.m[1][3] = 0.0; - nim->qto_xyz.m[2][0]=nim->qto_xyz.m[2][1]=nim->qto_xyz.m[2][3] = 0.0; - - /* last row is always [ 0 0 0 1 ] */ - - nim->qto_xyz.m[3][0]=nim->qto_xyz.m[3][1]=nim->qto_xyz.m[3][2] = 0.0; - nim->qto_xyz.m[3][3]= 1.0 ; - - nim->qform_code = NIFTI_XFORM_UNKNOWN ; - - } else { /** NIFTI: quaternion-specified transformation **/ - - nim->quatern_b = FIXED_FLOAT( nhdr.quatern_b ) ; - nim->quatern_c = FIXED_FLOAT( nhdr.quatern_c ) ; - nim->quatern_d = FIXED_FLOAT( nhdr.quatern_d ) ; - - nim->qoffset_x = FIXED_FLOAT(nhdr.qoffset_x) ; - nim->qoffset_y = FIXED_FLOAT(nhdr.qoffset_y) ; - nim->qoffset_z = FIXED_FLOAT(nhdr.qoffset_z) ; - - nim->qfac = (nhdr.pixdim[0] < 0.0) ? -1.0 : 1.0 ; /* left-handedness? */ - - nim->qto_xyz = quatern_to_mat44( - nim->quatern_b, nim->quatern_c, nim->quatern_c, - nim->qoffset_x, nim->qoffset_y, nim->qoffset_z, - nim->dx , nim->dy , nim->dz , - nim->qfac ) ; - - nim->qform_code = nhdr.qform_code ; - } - - /** load inverse transformation (x,y,z) -> (i,j,k) **/ - - nim->qto_ijk = mat44_inverse( nim->qto_xyz ) ; - - /** load sto_xyz affine transformation, if present **/ - - if( !is_nifti || nhdr.sform_code <= 0 ){ /** no sto transformation **/ - - nim->sform_code = NIFTI_XFORM_UNKNOWN ; - - } else { /** sto transformation from srow_*[] **/ - - nim->sto_xyz.m[0][0] = nhdr.srow_x[0] ; - nim->sto_xyz.m[0][1] = nhdr.srow_x[1] ; - nim->sto_xyz.m[0][2] = nhdr.srow_x[2] ; - nim->sto_xyz.m[0][3] = nhdr.srow_x[3] ; - - nim->sto_xyz.m[1][0] = nhdr.srow_y[0] ; - nim->sto_xyz.m[1][1] = nhdr.srow_y[1] ; - nim->sto_xyz.m[1][2] = nhdr.srow_y[2] ; - nim->sto_xyz.m[1][3] = nhdr.srow_y[3] ; - - nim->sto_xyz.m[2][0] = nhdr.srow_z[0] ; - nim->sto_xyz.m[2][1] = nhdr.srow_z[1] ; - nim->sto_xyz.m[2][2] = nhdr.srow_z[2] ; - nim->sto_xyz.m[2][3] = nhdr.srow_z[3] ; - - /* last row is always [ 0 0 0 1 ] */ - - nim->sto_xyz.m[3][0]=nim->sto_xyz.m[3][1]=nim->sto_xyz.m[3][2] = 0.0; - nim->sto_xyz.m[3][3]= 1.0 ; - - nim->sto_ijk = mat44_inverse( nim->sto_xyz ) ; - - nim->sform_code = nhdr.sform_code ; - } - - /* miscellaneous NIFTI stuff */ - - if( is_nifti ){ - nim->scl_slope = FIXED_FLOAT( nhdr.scl_slope ) ; - nim->scl_inter = FIXED_FLOAT( nhdr.scl_inter ) ; - - nim->intent_code = nhdr.intent_code ; - - nim->intent_p1 = FIXED_FLOAT( nhdr.intent_p1 ) ; - nim->intent_p2 = FIXED_FLOAT( nhdr.intent_p2 ) ; - nim->intent_p3 = FIXED_FLOAT( nhdr.intent_p3 ) ; - - nim->toffset = FIXED_FLOAT( nhdr.toffset ) ; - - memcpy(nim->intent_name,nhdr.intent_name,15); nim->intent_name[15] = '\0'; - - nim->xyz_units = XYZT_TO_SPACE(nhdr.xyzt_units) ; - nim->time_units = XYZT_TO_TIME (nhdr.xyzt_units) ; - - nim->freq_dim = DIM_INFO_TO_FREQ_DIM ( nhdr.dim_info ) ; - nim->phase_dim = DIM_INFO_TO_PHASE_DIM( nhdr.dim_info ) ; - nim->slice_dim = DIM_INFO_TO_SLICE_DIM( nhdr.dim_info ) ; - - nim->slice_code = nhdr.slice_code ; - nim->slice_start = nhdr.slice_start ; - nim->slice_end = nhdr.slice_end ; - nim->slice_duration = FIXED_FLOAT(nhdr.slice_duration) ; - } - - /* Miscellaneous ANALYZE stuff */ - - nim->cal_min = FIXED_FLOAT(nhdr.cal_min) ; - nim->cal_max = FIXED_FLOAT(nhdr.cal_max) ; - - memcpy(nim->descrip ,nhdr.descrip ,79) ; nim->descrip [79] = '\0' ; - memcpy(nim->aux_file,nhdr.aux_file,23) ; nim->aux_file[23] = '\0' ; - - /** read the data if desired, then bug out **/ - - nim->fname = strdup(hname) ; /* save input filename */ - nim->iname = iname ; /* save image filename */ + ioff = (int)nhdr.vox_offset ; + } nim->iname_offset = ioff ; - if( read_data ) nifti_image_load( nim ) ; + + /**- deal with file names if set */ + if (fname!=NULL) { + nifti_set_filenames(nim,fname,0,0); + if (nim->iname==NULL) { ERREX("bad filename"); } + } else { + nim->fname = NULL; + nim->iname = NULL; + } + + /* clear extension fields */ + nim->num_ext = 0; + nim->ext_list = NULL; + + return nim; +} + +#undef ERREX +#define ERREX(msg) \ + do{ fprintf(stderr,"** ERROR: nifti_image_open(%s): %s\n", \ + (hname != NULL) ? hname : "(null)" , (msg) ) ; \ + return fptr ; } while(0) + +/*************************************************************** + * nifti_image_open + ***************************************************************/ +/*! znzFile nifti_image_open( char *hname, char *opts , nifti_image **nim) + \brief Read in NIFTI-1 or ANALYZE-7.5 file (pair) header information into a nifti_image struct. + + - The image data is not read from disk (it may be read later using + nifti_image_load(), for example). + - The image data will be stored in whatever data format the + input data is; no scaling will be applied. + - DT_BINARY data is not supported. + - nifti_image_free() can be used to delete the returned struct, + when you are done with it. + + \param hname filename of dataset .hdr or .nii file + \param opts options string for opening the header file + \param nim pointer to pointer to nifti_image struct + (this routine allocates the nifti_image struct) + \return file pointer (gzippable) to the file with the image data, + ready for reading. + <br>NULL if something fails badly. + \sa nifti_image_load, nifti_image_free + */ +znzFile nifti_image_open(const char * hname, char * opts, nifti_image ** nim) +{ + znzFile fptr=NULL; + /* open the hdr and reading it in, but do not load the data */ + *nim = nifti_image_read(hname,0); + /* open the image file, ready for reading (compressed works for all reads) */ + if( ((*nim) == NULL) || ((*nim)->iname == NULL) || + ((*nim)->nbyper <= 0) || ((*nim)->nvox <= 0) ) + ERREX("bad header info") ; + + /* open image data file */ + fptr = znzopen( (*nim)->iname, opts, nifti_is_gzfile((*nim)->iname) ); + if( znz_isnull(fptr) ) ERREX("Can't open data file") ; + + return fptr; +} + + +/*----------------------------------------------------------------------*/ +/*! return an allocated and filled nifti_1_header struct + + Read the binary header from disk, and swap bytes if necessary. + + \return an allocated nifti_1_header struct, or NULL on failure + + \param hname name of file containing header + \param swapped if not NULL, return whether header bytes were swapped + \param check flag to check for invalid nifti_1_header + + \warning ASCII header type is not supported + + \sa nifti_image_read, nifti_image_free, nifti_image_read_bricks +*//*--------------------------------------------------------------------*/ +nifti_1_header * nifti_read_header(const char * hname, int * swapped, int check) +{ + nifti_1_header nhdr, * hptr; + znzFile fp; + int bytes, lswap; + char * hfile; + char fname[] = { "nifti_read_header" }; + + /* determine file name to use for header */ + hfile = nifti_findhdrname(hname); + if( hfile == NULL ){ + if( g_opts.debug > 0 ) + LNI_FERR(fname,"failed to find header file for", hname); + return NULL; + } else if( g_opts.debug > 1 ) + fprintf(stderr,"-d %s: found header filename '%s'\n",fname,hfile); + + fp = znzopen( hfile, "rb", nifti_is_gzfile(hfile) ); + if( znz_isnull(fp) ){ + if( g_opts.debug > 0 ) LNI_FERR(fname,"failed to open header file",hfile); + free(hfile); + return NULL; + } + + free(hfile); /* done with filename */ + + if( has_ascii_header(fp) == 1 ){ + znzclose( fp ); + if( g_opts.debug > 0 ) + LNI_FERR(fname,"ASCII header type not supported",hname); + return NULL; + } + + /* read the binary header */ + bytes = znzread( &nhdr, 1, sizeof(nhdr), fp ); + znzclose( fp ); /* we are done with the file now */ + + if( bytes < (int)sizeof(nhdr) ){ + if( g_opts.debug > 0 ){ + LNI_FERR(fname,"bad binary header read for file", hname); + fprintf(stderr," - read %d of %d bytes\n",bytes, (int)sizeof(nhdr)); + } + return NULL; + } + + /* now just decide on byte swapping */ + lswap = need_nhdr_swap(nhdr.dim[0], nhdr.sizeof_hdr); /* swap data flag */ + if( check && lswap < 0 ){ + LNI_FERR(fname,"bad nifti_1_header for file", hname); + return NULL; + } + + if( lswap ) { + if ( g_opts.debug > 3 ) disp_nifti_1_header("-d nhdr pre-swap: ", &nhdr); + swap_nifti_header( &nhdr , NIFTI_VERSION(nhdr) ) ; + } + + if ( g_opts.debug > 2 ) disp_nifti_1_header("-d nhdr post-swap: ", &nhdr); + + if ( check && ! nifti_hdr_looks_good(&nhdr) ){ + LNI_FERR(fname,"nifti_1_header looks bad for file", hname); + return NULL; + } + + /* all looks good, so allocate memory for and return the header */ + hptr = (nifti_1_header *)malloc(sizeof(nifti_1_header)); + if( ! hptr ){ + fprintf(stderr,"** nifti_read_hdr: failed to alloc nifti_1_header\n"); + return NULL; + } + + if( swapped ) *swapped = lswap; /* only if they care <sniff!> */ + + memcpy(hptr, &nhdr, sizeof(nifti_1_header)); + + return hptr; +} + + +/*----------------------------------------------------------------------*/ +/*! decide if this nifti_1_header structure looks reasonable + + Check dim[0], dim[1], sizeof_hdr, and datatype. + Check magic string for "n+1". + Maybe more tests will follow. + + \return 1 if the header seems valid, 0 otherwise + + \sa nifti_nim_is_valid, valid_nifti_extensions +*//*--------------------------------------------------------------------*/ +int nifti_hdr_looks_good(const nifti_1_header * hdr) +{ + int nbyper, swapsize; + int c, errs = 0; + + /* check dim[0] and sizeof_hdr */ + if( need_nhdr_swap(hdr->dim[0], hdr->sizeof_hdr) < 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** bad nhdr fields: dim0, sizeof_hdr = %d, %d\n", + hdr->dim[0], hdr->sizeof_hdr); + errs++; + } + + /* check the valid dimension sizes (maybe dim[0] is bad) */ + for( c = 1; c <= hdr->dim[0] && c <= 7; c++ ) + if( hdr->dim[c] <= 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** bad nhdr field: dim[%d] = %d\n",c,hdr->dim[c]); + errs++; + } + + /* check the magic string */ + if( (hdr->magic[0] != 'n') || + (hdr->magic[1] != 'i' && hdr->magic[1] != '+') || + (hdr->magic[2] != '1') || + (hdr->magic[3] != '\0') ) + { + if( g_opts.debug > 0 ) + fprintf(stderr, + "** bad nhdr field: magic = '%.4s', should be \"n+1\" or \"ni1\"\n" + " (in hex) magic = 0x%02x%02x%02x%02x\n" + " should be = 0x6e2b3100 or 0x6e693100\n", + hdr->magic, hdr->magic[0], hdr->magic[1], + hdr->magic[2], hdr->magic[3]); + errs++; + } + + /* check the datatype */ + if( hdr->datatype == DT_BINARY || hdr->datatype == DT_UNKNOWN ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** bad nhdr field: datatype = %d\n",hdr->datatype); + errs++; + } + + nifti_datatype_sizes(hdr->datatype, &nbyper, &swapsize); + if( nbyper == 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** bad nhdr field: datatype = %d\n",hdr->datatype); + errs++; + } + + if( errs ) return 0; /* problems */ + + if( g_opts.debug > 2 ) fprintf(stderr,"-d nifti header looks good\n"); + + return 1; /* looks good */ +} + + +/*---------------------------------------------------------------------- + * check whether byte swapping is needed + * + * dim[0] should be in [0,7], and sizeof_hdr should be accurate + * + * \returns > 0 : needs swap + * 0 : does not need swap + * < 0 : error condition + *----------------------------------------------------------------------*/ +static int need_nhdr_swap( short dim0, int hdrsize ) +{ + short d0 = dim0; /* so we won't have to swap them on the stack */ + int hsize = hdrsize; + + if( d0 != 0 ){ /* then use it for the check */ + if( d0 > 0 && d0 <= 7 ) return 0; + + nifti_swap_2bytes(1, &d0); /* swap? */ + if( d0 > 0 && d0 <= 7 ) return 1; + + if( g_opts.debug > 1 ){ + fprintf(stderr,"** bad swapped d0 = %d, unswapped = ", d0); + nifti_swap_2bytes(1, &d0); /* swap? */ + fprintf(stderr,"%d\n", d0); + } + + return -1; /* bad, naughty d0 */ + } + + /* dim[0] == 0 should not happen, but could, so try hdrsize */ + if( hsize == sizeof(nifti_1_header) ) return 0; + + nifti_swap_4bytes(1, &hsize); /* swap? */ + if( hsize == sizeof(nifti_1_header) ) return 1; + + if( g_opts.debug > 1 ){ + fprintf(stderr,"** bad swapped hsize = %d, unswapped = ", hsize); + nifti_swap_4bytes(1, &hsize); /* swap? */ + fprintf(stderr,"%d\n", hsize); + } + + return -2; /* bad, naughty hsize */ +} + + +/* use macro LNI_FILE_ERROR instead of ERREX() +#undef ERREX +#define ERREX(msg) \ + do{ fprintf(stderr,"** ERROR: nifti_image_read(%s): %s\n", \ + (hname != NULL) ? hname : "(null)" , (msg) ) ; \ + return NULL ; } while(0) +*/ + + +/*************************************************************** + * nifti_image_read + ***************************************************************/ +/*! \brief Read a nifti header and optionally the data, creating a nifti_image. + + - The data buffer will be byteswapped if necessary. + - The data buffer will not be scaled. + - The data buffer is allocated with calloc(). + + \param hname filename of the nifti dataset + \param read_data Flag, true=read data blob, false=don't read blob. + \return A pointer to the nifti_image data structure. + + \sa nifti_image_free, nifti_free_extensions, nifti_image_read_bricks +*/ +nifti_image *nifti_image_read( const char *hname , int read_data ) +{ + struct nifti_1_header nhdr ; + nifti_image *nim ; + znzFile fp ; + int rv, ii , filesize, remaining; + char fname[] = { "nifti_image_read" }; + char *hfile=NULL; + + if( g_opts.debug > 1 ){ + fprintf(stderr,"-d image_read from '%s', read_data = %d",hname,read_data); +#ifdef HAVE_ZLIB + fprintf(stderr,", HAVE_ZLIB = 1\n"); +#else + fprintf(stderr,", HAVE_ZLIB = 0\n"); +#endif + } + + /**- determine filename to use for header */ + hfile = nifti_findhdrname(hname); + if( hfile == NULL ){ + if(g_opts.debug > 0) + LNI_FERR(fname,"failed to find header file for", hname); + return NULL; /* check return */ + } else if( g_opts.debug > 1 ) + fprintf(stderr,"-d %s: found header filename '%s'\n",fname,hfile); + + if( nifti_is_gzfile(hfile) ) filesize = -1; /* unknown */ + else filesize = nifti_get_filesize(hfile); + + fp = znzopen(hfile, "rb", nifti_is_gzfile(hfile)); + if( znz_isnull(fp) ){ + if( g_opts.debug > 0 ) LNI_FERR(fname,"failed to open header file",hfile); + free(hfile); + return NULL; + } + + rv = has_ascii_header( fp ); + if( rv < 0 ){ + if( g_opts.debug > 0 ) LNI_FERR(fname,"short header read",hfile); + znzclose( fp ); + free(hfile); + return NULL; + } + else if ( rv == 1 ) /* process special file type */ + return nifti_read_ascii_image( fp, hfile, filesize, read_data ); + + /* else, just process normally */ + + /**- read binary header */ + + ii = znzread( &nhdr , 1 , sizeof(nhdr) , fp ) ; /* read the thing */ + + /* keep file open so we can check for exts. after nifti_convert_nhdr2nim() */ + + if( ii < sizeof(nhdr) ){ + if( g_opts.debug > 0 ){ + LNI_FERR(fname,"bad binary header read for file", hfile); + fprintf(stderr," - read %d of %d bytes\n",ii, (int)sizeof(nhdr)); + } + znzclose(fp) ; + free(hfile); + return NULL; + } + + /* create output image struct and set it up */ + + /**- convert all nhdr fields to nifti_image fields */ + nim = nifti_convert_nhdr2nim(nhdr,hfile); + + if( nim == NULL ){ + znzclose( fp ) ; /* close the file */ + if( g_opts.debug > 0 ) + LNI_FERR(fname,"cannot create nifti image from header",hfile); + free(hfile); /* had to save this for debug message */ + return NULL; + } + + if( g_opts.debug > 3 ){ + fprintf(stderr,"+d nifti_image_read(), have nifti image:\n"); + if( g_opts.debug > 2 ) nifti_image_infodump(nim); + } + + /**- check for extensions (any errors here means no extensions) */ + if( NIFTI_ONEFILE(nhdr) ) remaining = nim->iname_offset - sizeof(nhdr); + else remaining = filesize - sizeof(nhdr); + + (void)nifti_read_extensions(nim, fp, remaining); + + znzclose( fp ) ; /* close the file */ + free(hfile); + + /**- read the data if desired, then bug out */ + if( read_data ){ + if( nifti_image_load( nim ) < 0 ){ + nifti_image_free(nim); /* take ball, go home. */ + return NULL; + } + } + else nim->data = NULL ; + + return nim ; +} + + +/*---------------------------------------------------------------------- + * has_ascii_header - see if the NIFTI header is an ASCII format + * + * If the file starts with the ASCII string "<nifti_image", then + * process the dataset as a type-3 .nia file. + * + * return: -1 on error, 1 if true, or 0 if false + * + * NOTE: this is NOT part of the NIFTI-1 standard + *----------------------------------------------------------------------*/ +static int has_ascii_header( znzFile fp ) +{ + char buf[16]; + int nread; + + if( znz_isnull(fp) ) return 0; + + nread = znzread( buf, 1, 12, fp ); + buf[12] = '\0'; + + if( nread < 12 ) return -1; + + znzrewind(fp); /* move back to the beginning, and check */ + + if( strcmp(buf, "<nifti_image") == 0 ) return 1; + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! nifti_read_ascii_image - process as a type-3 .nia image file + + return NULL on failure + + NOTE: this is NOT part of the NIFTI-1 standard +*//*--------------------------------------------------------------------*/ +nifti_image * nifti_read_ascii_image(znzFile fp, char *fname, int flen, + int read_data) +{ + nifti_image * nim; + int slen, txt_size, remain, rv = 0; + char * sbuf, lfunc[25] = { "nifti_read_ascii_image" }; + + if( nifti_is_gzfile(fname) ){ + LNI_FERR(lfunc, "compressed file with negative offset", fname); + free(fname); znzclose(fp); return NULL; + } + slen = flen; /* slen will be our buffer length */ + + if( g_opts.debug > 1 ) + fprintf(stderr,"-d %s: have ASCII NIFTI file of size %d\n",fname,slen); + + if( slen > 65530 ) slen = 65530 ; + sbuf = (char *)calloc(sizeof(char),slen+1) ; + if( !sbuf ){ + fprintf(stderr,"** %s: failed to alloc %d bytes for sbuf",lfunc,65530); + free(fname); znzclose(fp); return NULL; + } + znzread( sbuf , 1 , slen , fp ) ; + nim = nifti_image_from_ascii( sbuf, &txt_size ) ; free( sbuf ) ; + if( nim == NULL ){ + LNI_FERR(lfunc,"failed nifti_image_from_ascii()",fname); + free(fname); znzclose(fp); return NULL; + } + nim->nifti_type = NIFTI_FTYPE_ASCII ; + + /* compute remaining space for extensions */ + remain = flen - txt_size - (int)nifti_get_volsize(nim); + if( remain > 4 ){ + /* read extensions (reposition file pointer, first) */ + znzseek(fp, txt_size, SEEK_SET); + (void) nifti_read_extensions(nim, fp, remain); + } + + free(fname); + znzclose( fp ) ; + + nim->iname_offset = -1 ; /* check from the end of the file */ + + if( read_data ) rv = nifti_image_load( nim ) ; else nim->data = NULL ; + /* check for nifti_image_load() failure, maybe bail out */ + if( read_data && rv != 0 ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d failed image_load, free nifti image struct\n"); + free(nim); + return NULL; + } + return nim ; } -/*--------------------------------------------------------------------------*/ -/* Load the image data from disk into an already-prepared image struct. -----------------------------------------------------------------------------*/ - -#undef ERREX -#define ERREX(msg) \ - do{ fprintf(stderr,"** ERROR: nifti_image_load: %s\n",(msg)) ; \ - return ; } while(0) - -void nifti_image_load( nifti_image *nim ) + +/*---------------------------------------------------------------------- + * Read the extensions into the nifti_image struct 08 Dec 2004 [rickr] + * + * This function is called just after the header struct is read in, and + * it is assumed the file pointer has not moved. The value in remain + * is assumed to be accurate, reflecting the bytes of space for potential + * extensions. + * + * return the number of extensions read in, or < 0 on error + *----------------------------------------------------------------------*/ +static int nifti_read_extensions( nifti_image *nim, znzFile fp, int remain ) +{ + nifti1_extender extdr; /* defines extension existence */ + nifti1_extension extn; /* single extension to process */ + nifti1_extension * Elist; /* list of processed extensions */ + int posn, count; + + if( !nim || znz_isnull(fp) ) { + if( g_opts.debug > 0 ) + fprintf(stderr,"** nifti_read_extensions: bad inputs (%p,%p)\n", + (void *)nim, (void *)fp); + return -1; + } + + posn = znztell(fp); + + if( (posn != sizeof(nifti_1_header)) && + (nim->nifti_type != NIFTI_FTYPE_ASCII) ) + fprintf(stderr,"** WARNING: posn not header size (%d, %d)\n", + posn, (int)sizeof(nifti_1_header)); + + if( g_opts.debug > 2 ) + fprintf(stderr,"-d nre: posn = %d, offset = %d, type = %d, remain = %d\n", + posn, nim->iname_offset, nim->nifti_type, remain); + + if( remain < 16 ){ + if( g_opts.debug > 2 ){ + if( g_opts.skip_blank_ext ) + fprintf(stderr,"-d no extender in '%s' is okay, as " + "skip_blank_ext is set\n",nim->fname); + else + fprintf(stderr,"-d no space for extensions\n"); + } + return 0; + } + + count = znzread( extdr.extension, 1, 4, fp ); /* get extender */ + + if( count < 4 ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d file '%s' is too short for an extender\n", + nim->fname); + return 0; + } + + if( extdr.extension[0] != 1 ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d extender[0] (%d) shows no extensions for '%s'\n", + extdr.extension[0], nim->fname); + return 0; + } + + remain -= 4; + if( g_opts.debug > 2 ) + fprintf(stderr,"-d found valid 4-byte extender, remain = %d\n", remain); + + /* so we expect extensions, but have no idea of how many there may be */ + + count = 0; + Elist = NULL; + while (nifti_read_next_extension(&extn, nim, remain, fp) > 0) + { + if( nifti_add_exten_to_list(&extn, &Elist, count+1) < 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** failed adding ext %d to list\n", count); + return -1; + } + + /* we have a new extension */ + if( g_opts.debug > 1 ){ + fprintf(stderr,"+d found extension #%d, code = 0x%x, size = %d\n", + count, extn.ecode, extn.esize); + if( extn.ecode == NIFTI_ECODE_AFNI && g_opts.debug > 2 ) /* ~XML */ + fprintf(stderr," AFNI extension: %.*s\n", + extn.esize-8,extn.edata); + else if( extn.ecode == NIFTI_ECODE_COMMENT && g_opts.debug > 2 ) + fprintf(stderr," COMMENT extension: %.*s\n", /* TEXT */ + extn.esize-8,extn.edata); + } + remain -= extn.esize; + count++; + } + + if( g_opts.debug > 2 ) fprintf(stderr,"+d found %d extension(s)\n", count); + + nim->num_ext = count; + nim->ext_list = Elist; + + return count; +} + + +/*----------------------------------------------------------------------*/ +/*! nifti_add_extension - add an extension, with a copy of the data + + Add an extension to the nim->ext_list array. + Fill this extension with a copy of the data, noting the + length and extension code. + + \param nim - nifti_image to add extension to + \param data - raw extension data + \param length - length of raw extension data + \param ecode - extension code + + \sa extension codes NIFTI_ECODE_* in nifti1_io.h + \sa nifti_free_extensions, valid_nifti_extensions, nifti_copy_extensions + + \return 0 on success, -1 on error (and free the entire list) +*//*--------------------------------------------------------------------*/ +int nifti_add_extension(nifti_image *nim, const char * data, int len, int ecode) +{ + nifti1_extension ext; + + /* error are printed in functions */ + if( nifti_fill_extension(&ext, data, len, ecode) ) return -1; + if( nifti_add_exten_to_list(&ext, &nim->ext_list, nim->num_ext+1)) return -1; + + nim->num_ext++; /* success, so increment */ + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/* nifti_add_exten_to_list - add a new nifti1_extension to the list + + We will append via "malloc, copy and free", because on an error, + the old data pointers must all be released (sorry realloc(), only + quality dolphins get to become part of St@rk!st brand tunafish). + + return 0 on success, -1 on error (and free the entire list) +*//*--------------------------------------------------------------------*/ +static int nifti_add_exten_to_list( nifti1_extension * new_ext, + nifti1_extension ** list, int new_length ) +{ + nifti1_extension * tmplist; + int count; + + tmplist = *list; + *list = (nifti1_extension *)malloc(new_length * sizeof(nifti1_extension)); + + /* check for failure first */ + if( ! *list ){ + fprintf(stderr,"** failed to alloc %d extension structs (%d bytes)\n", + new_length, new_length*(int)sizeof(nifti1_extension)); + if( !tmplist ) return -1; /* no old list to lose */ + + for ( count = 0; count < new_length-1; count++ ) + if( tmplist[count].edata ) free(tmplist[count].edata); + free(tmplist); + return -1; + } + + /* we have memory, so copy the old and insert the new */ + memcpy(*list, tmplist, (new_length-1)*sizeof(nifti1_extension)); + + /* for some reason, I just don't like struct copy... */ + (*list)[new_length-1].esize = new_ext->esize; + (*list)[new_length-1].ecode = new_ext->ecode; + (*list)[new_length-1].edata = new_ext->edata; + + if( g_opts.debug > 2 ) + fprintf(stderr,"+d allocated and appended extension #%d to list\n", + new_length); + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/* nifti_fill_extension - given data and length, fill an extension struct + + Allocate memory for data, copy data, set the size and code. + + return 0 on success, -1 on error (and free the entire list) +*//*--------------------------------------------------------------------*/ +static int nifti_fill_extension( nifti1_extension *ext, const char * data, + int len, int ecode) { - size_t ntot , ii , ioff ; - FILE *fp ; - + int esize; + + if( !ext || !data || len < 0 ){ + fprintf(stderr,"** fill_ext: bad params (%p,%p,%d)\n", + (void *)ext, data, len); + return -1; + } else if( ! nifti_is_valid_ecode(ecode) ){ + fprintf(stderr,"** fill_ext: invalid ecode %d\n", ecode); + return -1; + } + + /* compute esize, first : len+8, and take ceiling up to a mult of 16 */ + esize = len+8; + if( esize & 0xf ) esize = (esize + 0xf) & ~0xf; + ext->esize = esize; + + /* allocate esize-8 (maybe more than len), using calloc for fill */ + ext->edata = (char *)calloc(esize-8, sizeof(char)); + if( !ext->edata ){ + fprintf(stderr,"** NFE: failed to alloc %d bytes for extension\n",len); + return -1; + } + + memcpy(ext->edata, data, len); /* copy the data, using len */ + ext->ecode = ecode; /* set the ecode */ + + if( g_opts.debug > 2 ) + fprintf(stderr,"+d alloc %d bytes for ext len %d, ecode %d, esize %d\n", + esize-8, len, ecode, esize); + + return 0; +} + + +/*---------------------------------------------------------------------- + * nifti_read_next_extension - read a single extension from the file + * + * return (>= 0 is okay): + * + * success : esize + * no extension : 0 + * error : -1 + *----------------------------------------------------------------------*/ +static int nifti_read_next_extension( nifti1_extension * nex, nifti_image *nim, + int remain, znzFile fp ) +{ + int swap = nim->byteorder != nifti_short_order(); + int count, size, code; + + /* first clear nex */ + nex->esize = nex->ecode = 0; + nex->edata = NULL; + + if( remain < 16 ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d only %d bytes remain, so no extension\n", remain); + return 0; + } + + /* must start with 4-byte size and code */ + count = znzread( &size, 4, 1, fp ); + if( count == 1 ) count += znzread( &code, 4, 1, fp ); + + if( count != 2 ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d current extension read failed\n"); + znzseek(fp, -4*count, SEEK_CUR); /* back up past any read */ + return 0; /* no extension, no error condition */ + } + + if( swap ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d pre-swap exts: code %d, size %d\n", code, size); + + nifti_swap_4bytes(1, &size); + nifti_swap_4bytes(1, &code); + } + + if( g_opts.debug > 2 ) + fprintf(stderr,"-d potential extension: code %d, size %d\n", code, size); + + if( !nifti_check_extension(nim, size, code, remain) ){ + if( znzseek(fp, -8, SEEK_CUR) < 0 ){ /* back up past any read */ + fprintf(stderr,"** failure to back out of extension read!\n"); + return -1; + } + return 0; + } + + /* now get the actual data */ + nex->esize = size; + nex->ecode = code; + + size -= 8; /* subtract space for size and code in extension */ + nex->edata = (char *)malloc(size * sizeof(char)); + if( !nex->edata ){ + fprintf(stderr,"** failed to allocate %d bytes for extension\n",size); + return -1; + } + + count = znzread(nex->edata, 1, size, fp); + if( count < size ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"-d read only %d (of %d) bytes for extension\n", + count, size); + free(nex->edata); + nex->edata = NULL; + return -1; + } + + /* success! */ + if( g_opts.debug > 2 ) + fprintf(stderr,"+d successfully read extension, code %d, size %d\n", + nex->ecode, nex->esize); + + return nex->esize; +} + + +/*----------------------------------------------------------------------*/ +/*! for each extension, check code, size and data pointer +*//*--------------------------------------------------------------------*/ +int valid_nifti_extensions(const nifti_image * nim) +{ + nifti1_extension * ext; + int c, errs; + + if( nim->num_ext <= 0 || nim->ext_list == NULL ){ + if( g_opts.debug > 2 ) fprintf(stderr,"-d empty extension list\n"); + return 0; + } + + /* for each extension, check code, size and data pointer */ + ext = nim->ext_list; + errs = 0; + for ( c = 0; c < nim->num_ext; c++ ){ + if( ! nifti_is_valid_ecode(ext->ecode) ) { + if( g_opts.debug > 1 ) + fprintf(stderr,"-d ext %d, invalid code %d\n", c, ext->ecode); + errs++; + } + + if( ext->esize <= 0 ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d ext %d, bad size = %d\n", c, ext->esize); + errs++; + } else if( ext->esize & 0xf ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d ext %d, size %d not multiple of 16\n", + c, ext->esize); + errs++; + } + + if( ext->edata == NULL ){ + if( g_opts.debug > 1 ) fprintf(stderr,"-d ext %d, missing data\n", c); + errs++; + } + + ext++; + } + + if( errs > 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"-d had %d extension errors, none will be written\n", + errs); + return 0; + } + + /* if we're here, we're good */ + return 1; +} + + +/*----------------------------------------------------------------------*/ +/*! check whether the extension code is valid + + \return 1 if valid, 0 otherwise +*//*--------------------------------------------------------------------*/ +int nifti_is_valid_ecode( int ecode ) +{ + if( ecode < NIFTI_ECODE_IGNORE || /* minimum code number (0) */ + ecode > NIFTI_MAX_ECODE || /* maximum code number */ + ecode & 1 ) /* cannot be odd */ + return 0; + + return 1; +} + + +/*---------------------------------------------------------------------- + * check for valid size and code, as well as can be done + *----------------------------------------------------------------------*/ +static int nifti_check_extension(nifti_image *nim, int size, int code, int rem) +{ + /* check for bad code before bad size */ + if( ! nifti_is_valid_ecode(code) ) { + if( g_opts.debug > 2 ) + fprintf(stderr,"-d invalid extension code %d\n",code); + return 0; + } + + if( size < 16 ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d ext size %d, no extension\n",size); + return 0; + } + + if( size > rem ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d ext size %d, space %d, no extension\n", size, rem); + return 0; + } + + if( size & 0xf ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d nifti extension size %d not multiple of 16\n",size); + return 0; + } + + if( nim->nifti_type == NIFTI_FTYPE_ASCII && size > LNI_MAX_NIA_EXT_LEN ){ + if( g_opts.debug > 2 ) + fprintf(stderr,"-d NVE, bad nifti_type 3 size %d\n", size); + return 0; + } + + return 1; +} + + +/*---------------------------------------------------------------------- + * nifti_image_load_prep - prepare to read data + * + * Check nifti_image fields, open the file and seek to the appropriate + * offset for reading. + * + * return NULL on failure + *----------------------------------------------------------------------*/ +static znzFile nifti_image_load_prep( nifti_image *nim ) +{ + /* set up data space, open data file and seek, then call nifti_read_buffer */ + size_t ntot , ii , ioff; + znzFile fp; + char fname[] = { "nifti_image_load_prep" }; + + /**- perform sanity checks */ if( nim == NULL || nim->iname == NULL || - nim->nbyper <= 0 || nim->nvox <= 0 ) ERREX("bad input struct") ; - - /** open image data file **/ - - ntot = (size_t)(nim->nbyper) * (size_t)(nim->nvox) ; /* total bytes */ - - fp = fopen( nim->iname , "rb" ) ; - if( fp == NULL ) ERREX("Can't open data file") ; - if( nim->iname_offset < 0 ){ /* negative offset means */ - ii = get_filesize( nim->iname ) ; /* figure from end of file */ + nim->nbyper <= 0 || nim->nvox <= 0 ) + { + if ( g_opts.debug > 0 ){ + if( !nim ) fprintf(stderr,"** ERROR: N_image_load: no nifti image\n"); + else fprintf(stderr,"** ERROR: N_image_load: bad params (%p,%d,%d)\n", + nim->iname, nim->nbyper, nim->nvox); + } + return NULL; + } + + ntot = nifti_get_volsize(nim) ; /* total bytes to read */ + + /**- open image data file */ + + fp = znzopen(nim->iname, "rb", nifti_is_gzfile(nim->iname)); + if( znz_isnull(fp) ){ + if( g_opts.debug > 0 ) LNI_FERR(fname,"cannot open data file",nim->iname); + return NULL; + } + + /**- get image offset: a negative offset means to figure from end of file */ + if( nim->iname_offset < 0 ){ + if( nifti_is_gzfile(nim->iname) ){ + if( g_opts.debug > 0 ) + LNI_FERR(fname,"negative offset for compressed file",nim->iname); + znzclose(fp); + return NULL; + } + ii = nifti_get_filesize( nim->iname ) ; + if( ii <= 0 ){ + if( g_opts.debug > 0 ) LNI_FERR(fname,"empty data file",nim->iname); + znzclose(fp); + return NULL; + } ioff = (ii > ntot) ? ii-ntot : 0 ; } else { /* non-negative offset */ ioff = nim->iname_offset ; /* means use it directly */ } - fseek( fp , ioff , SEEK_SET ) ; - - /* make space for data, then read all of it in one operation */ - - if( nim->data != NULL ) free(nim->data) ; - - nim->data = (void *)calloc(1,ntot) ; - if( nim->data == NULL ) ERREX("can't malloc array space") ; - - ii = fread( nim->data , 1 , ntot , fp ) ; /*** data input! ***/ - fclose( fp ) ; - - /** if read was short, fill rest of array with 0 bytes **/ - + + /**- seek to the appropriate read position */ + if( znzseek(fp , ioff , SEEK_SET) < 0 ){ + fprintf(stderr,"** could not seek to offset %d in file '%s'\n", + (int)ioff, nim->iname); + znzclose(fp); + return NULL; + } + + /**- and return the File pointer */ + return fp; +} + + +/*---------------------------------------------------------------------- + * nifti_image_load + *----------------------------------------------------------------------*/ +/*! \fn int nifti_image_load( nifti_image *nim ) + \brief Load the image blob into a previously initialized nifti_image. + + - If not yet set, the data buffer is allocated with calloc(). + - The data buffer will be byteswapped if necessary. + - The data buffer will not be scaled. + + This function is used to read the image from disk. It should be used + after a function such as nifti_image_read(), so that the nifti_image + structure is already initialized. + + \param nim pointer to a nifti_image (previously initialized) + \return 0 on success, -1 on failure + \sa nifti_image_read, nifti_image_free, nifti_image_unload +*/ +int nifti_image_load( nifti_image *nim ) +{ + /* set up data space, open data file and seek, then call nifti_read_buffer */ + size_t ntot , ii ; + znzFile fp ; + + /**- open the file and position the FILE pointer */ + fp = nifti_image_load_prep( nim ); + + if( fp == NULL ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** nifti_image_load, failed load_prep\n"); + return -1; + } + + ntot = nifti_get_volsize(nim); + + /**- if the data pointer is not yet set, get memory space for the image */ + + if( nim->data == NULL ) + { + nim->data = (void *)calloc(1,ntot) ; /* create image memory */ + if( nim->data == NULL ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** failed to alloc %d bytes for image data\n", + (int)ntot); + znzclose(fp); + return -1; + } + } + + /**- now that everything is set up, do the reading */ + ii = nifti_read_buffer(fp,nim->data,ntot,nim); if( ii < ntot ){ - fprintf(stderr,"++ WARNING: nifti_image_load(%s):\n" - " data bytes needed = %u\n" - " data bytes input = %u\n" - " number missing = %u (set to 0)\n", - nim->iname , (unsigned int)ntot , - (unsigned int)ii , (unsigned int)(ntot-ii) ) ; - memset( (char *)(nim->data)+ii , 0 , ntot-ii ) ; - } - - /** byte swap array if needed **/ - - if( nim->swapsize > 1 && nim->byteorder != short_order() ) - swap_Nbytes( ntot / nim->swapsize , nim->swapsize , nim->data ) ; - + znzclose(fp) ; + free(nim->data) ; + nim->data = NULL ; + return -1 ; /* errors were printed in nifti_read_buffer() */ + } + + /**- close the file */ + znzclose( fp ) ; + + return 0 ; +} + + +/* 30 Nov 2004 [rickr] +#undef ERREX +#define ERREX(msg) \ + do{ fprintf(stderr,"** ERROR: nifti_read_buffer: %s\n",(msg)) ; \ + return 0; } while(0) +*/ + +/*----------------------------------------------------------------------*/ +/*! read ntot bytes of data from an open file and byte swaps if necessary + + note that nifti_image is required for information on datatype, bsize + (for any needed byte swapping), etc. + + This function does not allocate memory, so dataptr must be valid. +*//*--------------------------------------------------------------------*/ +size_t nifti_read_buffer(znzFile fp, void* dataptr, size_t ntot, + nifti_image *nim) +{ + size_t ii; + + if( dataptr == NULL ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** ERROR: nifti_read_buffer: NULL dataptr\n"); + return -1; + } + + ii = znzread( dataptr , 1 , ntot , fp ) ; /* data input */ + + /* if read was short, fail */ + if( ii < ntot ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"++ WARNING: nifti_read_buffer(%s):\n" + " data bytes needed = %u\n" + " data bytes input = %u\n" + " number missing = %u (set to 0)\n", + nim->iname , (unsigned int)ntot , + (unsigned int)ii , (unsigned int)(ntot-ii) ) ; + /* memset( (char *)(dataptr)+ii , 0 , ntot-ii ) ; now failure [rickr] */ + return -1 ; + } + + /* byte swap array if needed */ + + if( nim->swapsize > 1 && nim->byteorder != nifti_short_order() ) + nifti_swap_Nbytes( ntot / nim->swapsize , nim->swapsize , dataptr ) ; + #ifdef isfinite - /** check input float arrays for goodness, and fix bad numbers **/ - - switch( nim->datatype ){ - - case NIFTI_TYPE_FLOAT32: - case NIFTI_TYPE_COMPLEX64:{ - register float *far = (float *)nim->data ; register int jj,nj ; - nj = ntot / sizeof(float) ; - for( jj=0 ; jj < nj ; jj++ ) far[jj] = FIXED_FLOAT(far[jj]); - } - break ; - - case NIFTI_TYPE_FLOAT64: - case NIFTI_TYPE_COMPLEX128:{ - register double *far = (double *)nim->data ; register int jj,nj ; - nj = ntot / sizeof(double) ; - for( jj=0 ; jj < nj ; jj++ ) far[jj] = FIXED_FLOAT(far[jj]); - } - break ; - - } +{ + /* check input float arrays for goodness, and fix bad floats */ + int fix_count = 0 ; + + switch( nim->datatype ){ + + case NIFTI_TYPE_FLOAT32: + case NIFTI_TYPE_COMPLEX64:{ + register float *far = (float *)dataptr ; register int jj,nj ; + nj = ntot / sizeof(float) ; + for( jj=0 ; jj < nj ; jj++ ) /* count fixes 30 Nov 2004 [rickr] */ + if( !IS_GOOD_FLOAT(far[jj]) ){ + far[jj] = 0 ; + fix_count++ ; + } + } + break ; + + case NIFTI_TYPE_FLOAT64: + case NIFTI_TYPE_COMPLEX128:{ + register double *far = (double *)dataptr ; register int jj,nj ; + nj = ntot / sizeof(double) ; + for( jj=0 ; jj < nj ; jj++ ) /* count fixes 30 Nov 2004 [rickr] */ + if( !IS_GOOD_FLOAT(far[jj]) ){ + far[jj] = 0 ; + fix_count++ ; + } + } + break ; + + } + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d in image, %d bad floats were set to 0\n", fix_count); +} #endif - - return ; + + return ii; } /*--------------------------------------------------------------------------*/ -/* Unload the data in a nifti_image struct, but keep the metadata. -----------------------------------------------------------------------------*/ - +/*! Unload the data in a nifti_image struct, but keep the metadata. +*//*------------------------------------------------------------------------*/ void nifti_image_unload( nifti_image *nim ) { if( nim != NULL && nim->data != NULL ){ @@ -1478,102 +4399,275 @@ } /*--------------------------------------------------------------------------*/ -/* Free a nifti_image struct that was read by nifti_image_read(). -----------------------------------------------------------------------------*/ - +/*! free 'everything' about a nifti_image struct (including the passed struct) + + free (only fields which are not NULL): + - fname and iname + - data + - any ext_list[i].edata + - ext_list + - nim +*//*------------------------------------------------------------------------*/ void nifti_image_free( nifti_image *nim ) { if( nim == NULL ) return ; if( nim->fname != NULL ) free(nim->fname) ; if( nim->iname != NULL ) free(nim->iname) ; if( nim->data != NULL ) free(nim->data ) ; + (void)nifti_free_extensions( nim ) ; free(nim) ; return ; } + /*--------------------------------------------------------------------------*/ -/* Print to stdout some info about a nifti_image struct. -----------------------------------------------------------------------------*/ - -void nifti_image_infodump( nifti_image *nim ) +/*! free the nifti extensions + + - If any edata pointer is set in the extension list, free() it. + - Free ext_list, if it is set. + - Clear num_ext and ext_list from nim. + + \return 0 on success, -1 on error + + \sa nifti_add_extension, nifti_copy_extensions +*//*------------------------------------------------------------------------*/ +int nifti_free_extensions( nifti_image *nim ) { - char *str = nifti_image_to_ascii( nim ) ; - if( str != NULL ){ fputs(str,stdout) ; free(str) ; } - return ; + int c ; + if( nim == NULL ) return -1; + if( nim->num_ext > 0 && nim->ext_list ){ + for( c = 0; c < nim->num_ext; c++ ) + if ( nim->ext_list[c].edata ) free(nim->ext_list[c].edata); + free(nim->ext_list); + } + /* or if it is inconsistent, warn the user (if we are not in quiet mode) */ + else if ( (nim->num_ext > 0 || nim->ext_list != NULL) && (g_opts.debug > 0) ) + fprintf(stderr,"** warning: nifti extension num/ptr mismatch (%d,%p)\n", + nim->num_ext, (void *)nim->ext_list); + + if( g_opts.debug > 2 ) + fprintf(stderr,"+d free'd %d extension(s)\n", nim->num_ext); + + nim->num_ext = 0; + nim->ext_list = NULL; + + return 0; } + /*--------------------------------------------------------------------------*/ -/* Write a nifti_image to disk. The following fields of nim affect how - the output appears: - - nifti_type = 0 ==> ANALYZE-7.5 format file pair will be written - - nifti_type = 1 ==> NIFTI-1 format single file will be written - (data offset will be 348) - - nifti_type = 2 ==> NIFTI_1 format file pair will be written - - nifti_type = 3 ==> NIFTI_1 ASCII single file will be written - - fname is the name of the output file (header or header+data) - - if a file pair is being written, iname is the name of the data file - - existing files WILL be overwritten with extreme prejudice - - if qform_code > 0, the quatern_*, qoffset_*, and qfac fields determine - the qform output, NOT the qto_xyz matrix; if you want to compute these - fields from the qto_xyz matrix, you can use the utility function - mat44_to_quatern() -----------------------------------------------------------------------------*/ - -#undef ERREX -#define ERREX(msg) \ - do{ fprintf(stderr,"** ERROR: nifti_image_write: %s\n",(msg)) ; \ - return ; } while(0) - -void nifti_image_write( nifti_image *nim ) +/*! Print to stdout some info about a nifti_image struct. +*//*------------------------------------------------------------------------*/ +void nifti_image_infodump( const nifti_image *nim ) +{ + char *str = nifti_image_to_ascii( nim ) ; + /* stdout -> stderr 2 Dec 2004 [rickr] */ + if( str != NULL ){ fputs(str,stderr) ; free(str) ; } + return ; +} + + +/*-------------------------------------------------------------------------- + * nifti_write_buffer just check for a null znzFile and call znzwrite + *--------------------------------------------------------------------------*/ +/*! \fn size_t nifti_write_buffer(znzFile fp, void *buffer, size_t numbytes) + \brief write numbytes of buffer to file, fp + + \param fp File pointer (from znzopen) to gzippable nifti datafile + \param buffer data buffer to be written + \param numbytes number of bytes in buffer to write + \return number of bytes successfully written +*/ +size_t nifti_write_buffer(znzFile fp, const void *buffer, size_t numbytes) +{ + /* Write all the image data at once (no swapping here) */ + size_t ss; + if (znz_isnull(fp)){ + fprintf(stderr,"** ERROR: nifti_write_buffer: null file pointer\n"); + return 0; + } + ss = znzwrite( (void*)buffer , 1 , numbytes , fp ) ; + return ss; +} + + +/*----------------------------------------------------------------------*/ +/*! write the nifti_image data to file (from nim->data or from NBL) + + If NBL is not NULL, write the data from that structure. Otherwise, + write it out from nim->data. No swapping is done here. + + \param fp : File pointer + \param nim : nifti_image corresponding to the data + \param NBL : optional source of write data (if NULL use nim->data) + + \return 0 on success, -1 on failure + + Note: the nifti_image byte_order is set as that of the current CPU. + This is because such a conversion was made to the data upon + reading, while byte_order was not set (so the programs would + know what format the data was on disk). Effectively, since + byte_order should match what is on disk, it should bet set to + that of the current CPU whenever new filenames are assigned. +*//*--------------------------------------------------------------------*/ +int nifti_write_all_data(znzFile fp, nifti_image * nim, + const nifti_brick_list * NBL) { - struct nifti_1_header nhdr ; - FILE *fp ; - size_t ss ; - - if( nim == NULL ) ERREX("NULL input") ; - if( nim->fname == NULL || nim->fname[0] == '\0' ) ERREX("bad fname input") ; - if( nim->data == NULL ) ERREX("no image data") ; - - /* make iname from fname, if needed */ - - switch( nim->nifti_type ){ - - default: /* writing into 2 files */ - if( nim->iname != NULL && strcmp(nim->fname,nim->fname) == 0 ){ - free(nim->iname) ; nim->iname = NULL ; - } - if( nim->iname == NULL ){ - int ll = strlen(nim->fname) ; - nim->iname = (char *)calloc(1,ll+5) ; - strcpy(nim->iname,nim->fname) ; - if( ll > 4 ) strcpy(nim->iname+ll-4,".img") ; /* create .img filename */ - else strcat(nim->iname ,".img") ; - } - nim->iname_offset = 0 ; - break ; - - case 1: /* NIFTI-1 single binary file */ - nim->iname_offset = sizeof(nhdr) ; - break ; - - /* non-standard case: */ - case 3:{ /* NIFTI-1 ASCII header + binary data (single file) */ - char *hstr ; - nim->iname_offset = -1 ; /* compute offset from filesize */ - nim->byteorder = short_order() ; /* am writing in current order */ - hstr = nifti_image_to_ascii( nim ) ; /* get header in ASCII form */ - if( hstr == NULL ) ERREX("bad ASCII header creation?") ; - fp = fopen( nim->fname , "wb" ) ; - if( fp == NULL ){ free(hstr); ERREX("can't open output file"); } - fputs(hstr,fp) ; - } - goto DataWriter ; /* writes the binary data to fp */ - } - - /***** Here -- write a binary header *****/ + size_t ss; + int bnum; + + if( !NBL ){ /* just write one buffer and get out of here */ + if( nim->data == NULL ){ + fprintf(stderr,"** NWAD: no image data to write\n"); + return -1; + } + + ss = nifti_write_buffer(fp,nim->data,nim->nbyper * nim->nvox); + if (ss < (size_t)(nim->nbyper * nim->nvox)){ + fprintf(stderr, + "** ERROR: NWAD: wrote only %d of %d bytes to file\n", + (int)ss, nim->nbyper * nim->nvox); + return -1; + } + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d wrote single image of %d bytes\n",(int)ss); + } else { + if( ! NBL->bricks || NBL->nbricks <= 0 || NBL->bsize <= 0 ){ + fprintf(stderr,"** NWAD: no brick data to write (%p,%d,%d)\n", + (void *)NBL->bricks, NBL->nbricks, NBL->bsize); + return -1; + } + + for( bnum = 0; bnum < NBL->nbricks; bnum++ ){ + ss = nifti_write_buffer(fp, NBL->bricks[bnum], NBL->bsize); + if( ss < (size_t)NBL->bsize ){ + fprintf(stderr, + "** NWAD ERROR: wrote %d of %d bytes of brick %d of %d to file", + (int)ss, NBL->bsize, bnum+1, NBL->nbricks); + return -1; + } + } + if( g_opts.debug > 1 ) + fprintf(stderr,"+d wrote image of %d brick(s), each of %d bytes\n", + NBL->nbricks, NBL->bsize); + } + + /* mark as being in this CPU byte order */ + nim->byteorder = nifti_short_order() ; + + return 0; +} + +/* return number of extensions written, or -1 on error */ +static int nifti_write_extensions(znzFile fp, nifti_image *nim) +{ + nifti1_extension * list; + char extdr[4] = { 0, 0, 0, 0 }; + int c, size, ok = 1; + + if( znz_isnull(fp) || !nim || nim->num_ext < 0 ){ + if( g_opts.debug > 0 ) + fprintf(stderr,"** nifti_write_extensions, bad params\n"); + return -1; + } + + /* if no extensions and user requests it, skip extender */ + if( g_opts.skip_blank_ext && (nim->num_ext == 0 || ! nim->ext_list ) ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"-d no exts and skip_blank_ext set, " + "so skipping 4-byte extender\n"); + return 0; + } + + /* if invalid extension list, clear num_ext */ + if( ! valid_nifti_extensions(nim) ) nim->num_ext = 0; + + /* write out extender block */ + if( nim->num_ext > 0 ) extdr[0] = 1; + if( nifti_write_buffer(fp, extdr, 4) != 4 ){ + fprintf(stderr,"** failed to write extender\n"); + return -1; + } + + list = nim->ext_list; + for ( c = 0; c < nim->num_ext; c++ ){ + size = nifti_write_buffer(fp, &list->esize, sizeof(int)); + ok = (size == (int)sizeof(int)); + if( ok ){ + size = nifti_write_buffer(fp, &list->ecode, sizeof(int)); + ok = (size == (int)sizeof(int)); + } + if( ok ){ + size = nifti_write_buffer(fp, list->edata, list->esize - 8); + ok = (size == list->esize - 8); + } + + if( !ok ){ + fprintf(stderr,"** failed while writing extension #%d\n",c); + return -1; + } + + list++; + } + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d wrote out %d extension(s)\n", nim->num_ext); + + return nim->num_ext; +} + + +/*----------------------------------------------------------------------*/ +/*! basic initialization of a nifti_image struct (to a 1x1x1 image) +*//*--------------------------------------------------------------------*/ +nifti_image* nifti_simple_init_nim(void) +{ + nifti_image *nim; + struct nifti_1_header nhdr; + int nbyper, swapsize; + + memset(&nhdr,0,sizeof(nhdr)) ; /* zero out header, to be safe */ + + nhdr.sizeof_hdr = sizeof(nhdr) ; + nhdr.regular = 'r' ; /* for some stupid reason */ + + nhdr.dim[0] = 3 ; + nhdr.dim[1] = 1 ; nhdr.dim[2] = 1 ; nhdr.dim[3] = 1 ; + nhdr.dim[4] = 0 ; + + nhdr.pixdim[0] = 0.0 ; + nhdr.pixdim[1] = 1.0 ; nhdr.pixdim[2] = 1.0 ; + nhdr.pixdim[3] = 1.0 ; + + nhdr.datatype = DT_FLOAT32 ; + nifti_datatype_sizes( nhdr.datatype , &nbyper, &swapsize ); + nhdr.bitpix = 8 * nbyper ; + + nim = nifti_convert_nhdr2nim(nhdr,NULL); + nim->fname = NULL; + nim->iname = NULL; + return nim; +} + + +/*----------------------------------------------------------------------*/ +/*! convert a nifti_image structure to a nifti_1_header struct + + No allocation is done, this should be used via structure copy. + As in: + <pre> + nifti_1_header my_header; + my_header = nifti_convert_nim2nhdr(my_nim_pointer); + </pre> +*//*--------------------------------------------------------------------*/ +struct nifti_1_header nifti_convert_nim2nhdr(const nifti_image * nim) +{ + struct nifti_1_header nhdr; memset(&nhdr,0,sizeof(nhdr)) ; /* zero out header, to be safe */ - /** load the ANALYZE-7.5 generic parts of the header struct **/ + + /**- load the ANALYZE-7.5 generic parts of the header struct */ nhdr.sizeof_hdr = sizeof(nhdr) ; nhdr.regular = 'r' ; /* for some stupid reason */ @@ -1609,12 +4703,17 @@ memcpy(nhdr.aux_file ,nim->aux_file ,23) ; nhdr.aux_file[23] = '\0' ; } - /** Load NIFTI specific stuff into the header **/ - - if( nim->nifti_type > 0 ){ - - if( nim->nifti_type == 1 ) strcpy(nhdr.magic,"n+1") ; /* 1 file */ - else strcpy(nhdr.magic,"ni1") ; /* 2 files */ + /**- Load NIFTI specific stuff into the header */ + + if( nim->nifti_type > NIFTI_FTYPE_ANALYZE ){ /* then not ANALYZE */ + + if( nim->nifti_type == NIFTI_FTYPE_NIFTI1_1 ) strcpy(nhdr.magic,"n+1") ; + else strcpy(nhdr.magic,"ni1") ; + + nhdr.pixdim[1] = fabs(nhdr.pixdim[1]) ; nhdr.pixdim[2] = fabs(nhdr.pixdim[2]) ; + nhdr.pixdim[3] = fabs(nhdr.pixdim[3]) ; nhdr.pixdim[4] = fabs(nhdr.pixdim[4]) ; + nhdr.pixdim[5] = fabs(nhdr.pixdim[5]) ; nhdr.pixdim[6] = fabs(nhdr.pixdim[6]) ; + nhdr.pixdim[7] = fabs(nhdr.pixdim[7]) ; nhdr.intent_code = nim->intent_code ; nhdr.intent_p1 = nim->intent_p1 ; @@ -1664,35 +4763,396 @@ nhdr.slice_duration = nim->slice_duration ; } - /** Open file, write header **/ - - fp = fopen( nim->fname , "wb" ) ; - if( fp == NULL ) ERREX("can't open output file") ; - - ss = fwrite( &nhdr , 1 , sizeof(nhdr) , fp ) ; + return nhdr; +} + + +/*----------------------------------------------------------------------*/ +/*! \fn int nifti_copy_extensions(nifti_image * nim_dest, nifti_image * nim_src) + \brief copy the nifti1_extension list from src to dest + + Duplicate the list of nifti1_extensions. The dest structure must + be clear of extensions. + \return 0 on success, -1 on failure + + \sa nifti_add_extension, nifti_free_extensions +*/ +int nifti_copy_extensions(nifti_image * nim_dest, const nifti_image * nim_src) +{ + char * data; + size_t bytes; + int c, size, old_size; + + if( nim_dest->num_ext > 0 || nim_dest->ext_list != NULL ){ + fprintf(stderr,"** will not copy extensions over existing ones\n"); + return -1; + } + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d duplicating %d extension(s)\n", nim_src->num_ext); + + if( nim_src->num_ext <= 0 ) return 0; + + bytes = nim_src->num_ext * sizeof(nifti1_extension); /* I'm lazy */ + nim_dest->ext_list = (nifti1_extension *)malloc(bytes); + if( !nim_dest->ext_list ){ + fprintf(stderr,"** failed to allocate %d nifti1_extension structs\n", + nim_src->num_ext); + return -1; + } + + /* copy the extension data */ + nim_dest->num_ext = 0; + for( c = 0; c < nim_src->num_ext; c++ ){ + size = old_size = nim_src->ext_list[c].esize; + if( size & 0xf ) size = (size + 0xf) & ~0xf; /* make multiple of 16 */ + if( g_opts.debug > 2 ) + fprintf(stderr,"+d dup'ing ext #%d of size %d (from size %d)\n", + c, size, old_size); + data = (char *)calloc(size,sizeof(char)); /* calloc, maybe size > old */ + if( !data ){ + fprintf(stderr,"** failed to alloc %d bytes for extention\n", size); + if( c == 0 ) { free(nim_dest->ext_list); nim_dest->ext_list = NULL; } + /* otherwise, keep what we have (a.o.t. deleting them all) */ + return -1; + } + /* finally, fill the new structure */ + nim_dest->ext_list[c].esize = size; + nim_dest->ext_list[c].ecode = nim_src->ext_list[c].ecode; + nim_dest->ext_list[c].edata = data; + memcpy(data, nim_src->ext_list[c].edata, old_size); + + nim_dest->num_ext++; + } + + return 0; +} + + +/*----------------------------------------------------------------------*/ +/*! compute the total size of all extensions + + \return the total of all esize fields + + Note that each esize includes 4 bytes for ecode, 4 bytes for esize, + and the bytes used for the data. Each esize also needs to be a + multiple of 16, so it may be greater than the sum of its 3 parts. +*//*--------------------------------------------------------------------*/ +int nifti_extension_size(nifti_image *nim) +{ + int c, size = 0; + + if( !nim || nim->num_ext <= 0 ) return 0; + + if( g_opts.debug > 2 ) fprintf(stderr,"-d ext sizes:"); + + for ( c = 0; c < nim->num_ext; c++ ){ + size += nim->ext_list[c].esize; + if( g_opts.debug > 2 ) fprintf(stderr," %d",nim->ext_list[c].esize); + } + + if( g_opts.debug > 2 ) fprintf(stderr," (total = %d)\n",size); + + return size; +} + + +/*----------------------------------------------------------------------*/ +/*! set the nifti_image iname_offset field, based on nifti_type + + - if writing to 2 files, set offset to 0 + - if writing to a single NIFTI-1 file, set the offset to + 352 + total extension size, then align to 16-byte boundary + - if writing an ASCII header, set offset to -1 +*//*--------------------------------------------------------------------*/ +void nifti_set_iname_offset(nifti_image *nim) +{ + int offset; + + switch( nim->nifti_type ){ + + default: /* writing into 2 files */ + /* we only write files with 0 offset in the 2 file format */ + nim->iname_offset = 0 ; + break ; + + /* NIFTI-1 single binary file - always update */ + case NIFTI_FTYPE_NIFTI1_1: + offset = nifti_extension_size(nim)+sizeof(struct nifti_1_header)+4; + /* be sure offset is aligned to a 16 byte boundary */ + if ( ( offset % 16 ) != 0 ) offset = ((offset + 0xf) & ~0xf); + if( nim->iname_offset != offset ){ + if( g_opts.debug > 1 ) + fprintf(stderr,"+d changing offset from %d to %d\n", + nim->iname_offset, offset); + nim->iname_offset = offset; + } + break ; + + /* non-standard case: NIFTI-1 ASCII header + binary data (single file) */ + case NIFTI_FTYPE_ASCII: + nim->iname_offset = -1 ; /* compute offset from filesize */ + break ; + } +} + + +/*----------------------------------------------------------------------*/ +/*! write the nifti_image dataset to disk, optionally including data + + This is just a front-end for nifti_image_write_hdr_img2. + + \param nim nifti_image to write to disk + \param write_data write options (see nifti_image_write_hdr_img2) + \param opts file open options ("wb" from nifti_image_write) + + \sa nifti_image_write, nifti_image_write_hdr_img2, nifti_image_free, + nifti_set_filenames +*//*--------------------------------------------------------------------*/ +znzFile nifti_image_write_hdr_img( nifti_image *nim , int write_data , + const char* opts ) +{ + return nifti_image_write_hdr_img2(nim,write_data,opts,NULL,NULL); +} + + +#undef ERREX +#define ERREX(msg) \ + do{ fprintf(stderr,"** ERROR: nifti_image_write_hdr_img: %s\n",(msg)) ; \ + return fp ; } while(0) + + +/* ----------------------------------------------------------------------*/ +/*! This writes the header (and optionally the image data) to file + * + * If the image data file is left open it returns a valid znzFile handle. + * It also uses imgfile as the open image file is not null, and modifies + * it inside. + * + * \param nim nifti_image to write to disk + * \param write_opts flags whether to write data and/or close file (see below) + * \param opts file-open options, probably "wb" from nifti_image_write() + * \param imgfile optional open znzFile struct, for writing image data + (may be NULL) + * \param NBL optional nifti_brick_list, containing the image data + (may be NULL) + * + * Values for write_opts mode are based on two binary flags + * ( 0/1 for no-write/write data, and 0/2 for close/leave-open files ) : + * - 0 = do not write data and close (do not open data file) + * - 1 = write data and close + * - 2 = do not write data and leave data file open + * - 3 = write data and leave data file open + * + * \sa nifti_image_write, nifti_image_write_hdr_img, nifti_image_free, + * nifti_set_filenames +*//*---------------------------------------------------------------------*/ +znzFile nifti_image_write_hdr_img2(nifti_image *nim, int write_opts, + const char * opts, znzFile imgfile, const nifti_brick_list * NBL) +{ + struct nifti_1_header nhdr ; + znzFile fp=NULL; + size_t ss ; + int write_data, leave_open; + char func[] = { "nifti_image_write_hdr_img2" }; + + write_data = write_opts & 1; /* just separate the bits now */ + leave_open = write_opts & 2; + + if( ! nim ) ERREX("NULL input") ; + if( ! nifti_validfilename(nim->fname) ) ERREX("bad fname input") ; + if( write_data && ! nim->data && ! NBL ) ERREX("no image data") ; + + nifti_set_iname_offset(nim); + + if( g_opts.debug > 1 ){ + fprintf(stderr,"-d writing nifti file '%s'...\n", nim->fname); + if( g_opts.debug > 2 ) + fprintf(stderr,"-d nifti type %d, offset %d\n", + nim->nifti_type, nim->iname_offset); + } + + if( nim->nifti_type == NIFTI_FTYPE_ASCII ) /* non-standard case */ + return nifti_write_ascii_image(nim,NBL,opts,write_data,leave_open); + + nhdr = nifti_convert_nim2nhdr(nim); /* create the nifti1_header struct */ + + /* if writing to 2 files, make sure iname is set and different from fname */ + if( nim->nifti_type != NIFTI_FTYPE_NIFTI1_1 ){ + if( nim->iname && strcmp(nim->iname,nim->fname) == 0 ){ + free(nim->iname) ; nim->iname = NULL ; + } + if( nim->iname == NULL ){ /* then make a new one */ + nim->iname = nifti_makeimgname(nim->fname,nim->nifti_type,0,0); + if( nim->iname == NULL ) return NULL; + } + } + + /* if we have an imgfile and will write the header there, use it */ + if( ! znz_isnull(imgfile) && nim->nifti_type == NIFTI_FTYPE_NIFTI1_1 ){ + if( g_opts.debug > 2 ) fprintf(stderr,"+d using passed file for hdr\n"); + fp = imgfile; + } + else { + if( g_opts.debug > 2 ) + fprintf(stderr,"+d opening output file '%s'\n",nim->fname); + fp = znzopen( nim->fname , opts , nifti_is_gzfile(nim->fname) ) ; + if( znz_isnull(fp) ){ + LNI_FERR(func,"cannot open output file",nim->fname); + return fp; + } + } + + /* write the header and extensions */ + + ss = znzwrite(&nhdr , 1 , sizeof(nhdr) , fp); /* write header */ if( ss < sizeof(nhdr) ){ - fclose(fp) ; ERREX("bad write to output file") ; - } - - /** If not writing 1 file, close header and open image file **/ - - if( nim->nifti_type != 1 ){ - fclose(fp) ; - fp = fopen( nim->iname , "wb" ) ; - if( fp == NULL ) ERREX("can't open image file") ; - } - - /** Write all the image data at once (no swapping here) **/ - - DataWriter: - ss = fwrite( nim->data , nim->nbyper , nim->nvox , fp ) ; - fclose(fp) ; - if( ss < nim->nvox ) ERREX("incomplete write to image file") ; - - nim->byteorder = short_order() ; /* mark as being in this CPU byte order */ - return ; + LNI_FERR(func,"bad header write to output file",nim->fname); + znzclose(fp); return fp; + } + + /* partial file exists, and errors have been printed, so ignore return */ + if( nim->nifti_type != NIFTI_FTYPE_ANALYZE ) + (void)nifti_write_extensions(fp,nim); + + /* if the header is all we want, we are done */ + if( ! write_data && ! leave_open ){ + if( g_opts.debug > 2 ) fprintf(stderr,"-d header is all we want: done\n"); + znzclose(fp); return(fp); + } + + if( nim->nifti_type != NIFTI_FTYPE_NIFTI1_1 ){ /* get a new file pointer */ + znzclose(fp); /* first, close header file */ + if( ! znz_isnull(imgfile) ){ + if(g_opts.debug > 2) fprintf(stderr,"+d using passed file for img\n"); + fp = imgfile; + } + else { + if( g_opts.debug > 2 ) + fprintf(stderr,"+d opening img file '%s'\n", nim->iname); + fp = znzopen( nim->iname , opts , nifti_is_gzfile(nim->iname) ) ; + if( znz_isnull(fp) ) ERREX("cannot open image file") ; + } + } + + znzseek(fp, nim->iname_offset, SEEK_SET); /* in any case, seek to offset */ + + if( write_data ) nifti_write_all_data(fp,nim,NBL); + if( ! leave_open ) znzclose(fp); + + return fp; +} + + +/*----------------------------------------------------------------------*/ +/*! write a nifti_image to disk in ASCII format +*//*--------------------------------------------------------------------*/ +znzFile nifti_write_ascii_image(nifti_image *nim, const nifti_brick_list * NBL, + const char *opts, int write_data, int leave_open) +{ + znzFile fp; + char * hstr; + + hstr = nifti_image_to_ascii( nim ) ; /* get header in ASCII form */ + if( ! hstr ){ fprintf(stderr,"** failed image_to_ascii()\n"); return NULL; } + + fp = znzopen( nim->fname , opts , nifti_is_gzfile(nim->fname) ) ; + if( znz_isnull(fp) ){ + free(hstr); + fprintf(stderr,"** failed to open '%s' for ascii write\n",nim->fname); + return fp; + } + + znzputs(hstr,fp); /* header */ + nifti_write_extensions(fp,nim); /* extensions */ + + if ( write_data ) { nifti_write_all_data(fp,nim,NBL); } /* data */ + if ( ! leave_open ) { znzclose(fp); } + free(hstr); + return fp; /* returned but may be closed */ } + +/*--------------------------------------------------------------------------*/ +/*! Write a nifti_image to disk. + + Since data is properly byte-swapped upon reading, it is assumed + to be in the byte-order of the current CPU at write time. Thus, + nim->byte_order should match that of the current CPU. Note that + the nifti_set_filenames() function takes the flag, set_byte_order. + + The following fields of nim affect how the output appears: + - nifti_type = 0 ==> ANALYZE-7.5 format file pair will be written + - nifti_type = 1 ==> NIFTI-1 format single file will be written + (data offset will be 352+extensions) + - nifti_type = 2 ==> NIFTI_1 format file pair will be written + - nifti_type = 3 ==> NIFTI_1 ASCII single file will be written + - fname is the name of the output file (header or header+data) + - if a file pair is being written, iname is the name of the data file + - existing files WILL be overwritten with extreme prejudice + - if qform_code > 0, the quatern_*, qoffset_*, and qfac fields determine + the qform output, NOT the qto_xyz matrix; if you want to compute these + fields from the qto_xyz matrix, you can use the utility function + nifti_mat44_to_quatern() + + \sa nifti_image_write_bricks, nifti_image_free, nifti_set_filenames, + nifti_image_write_hdr_img +*//*------------------------------------------------------------------------*/ +void nifti_image_write( nifti_image *nim ) +{ + znzFile fp = nifti_image_write_hdr_img(nim,1,"wb"); + if( fp ){ + if( g_opts.debug > 2 ) fprintf(stderr,"-d niw: done with znzFile\n"); + free(fp); + } + if( g_opts.debug > 1 ) fprintf(stderr,"-d nifti_image_write: done\n"); +} + + +/*----------------------------------------------------------------------*/ +/*! similar to nifti_image_write, but data is in NBL struct, not nim->data + + \sa nifti_image_write, nifti_image_free, nifti_set_filenames, nifti_free_NBL +*//*--------------------------------------------------------------------*/ +void nifti_image_write_bricks( nifti_image *nim, const nifti_brick_list * NBL ) +{ + znzFile fp = nifti_image_write_hdr_img2(nim,1,"wb",NULL,NBL); + if( fp ){ + if( g_opts.debug > 2 ) fprintf(stderr,"-d niwb: done with znzFile\n"); + free(fp); + } + if( g_opts.debug > 1 ) fprintf(stderr,"-d niwb: done writing bricks\n"); +} + + +/*----------------------------------------------------------------------*/ +/*! copy the nifti_image structure, without data + + Duplicate the structure, including fname, iname and extensions. + Leave the data pointer as NULL. +*//*--------------------------------------------------------------------*/ +nifti_image * nifti_copy_nim_info(const nifti_image * src) +{ + nifti_image *dest; + dest = (nifti_image *)calloc(1,sizeof(nifti_image)); + if( !dest ){ + fprintf(stderr,"** NCNI: failed to alloc nifti_image\n"); + return NULL; + } + memcpy(dest, src, sizeof(nifti_image)); + if( src->fname ) dest->fname = nifti_strdup(src->fname); + if( src->iname ) dest->iname = nifti_strdup(src->iname); + dest->num_ext = 0; + dest->ext_list = NULL; + /* errors will be printed in NCE(), continue in either case */ + (void)nifti_copy_extensions(dest, src); + + dest->data = NULL; + + return dest; +} + + /*------------------------------------------------------------------------*/ /* Un-escape a C string in place -- that is, convert XML escape sequences back into their characters. (This can be done in place since the @@ -1712,7 +5172,7 @@ #define CR 0x0D #define LF 0x0A -int unescape_string( char *str ) +static int unescape_string( char *str ) { int ii,jj , nn,ll ; @@ -1809,13 +5269,13 @@ The result should be free()-ed when you are done with it. --------------------------------------------------------------------------*/ -char *escapize_string( char *str ) +static char *escapize_string( const char * str ) { int ii,jj , lstr,lout ; char *out ; if( str == NULL || (lstr=strlen(str)) == 0 ){ /* 0 length */ - out = strdup("''") ; return out ; /* string?? */ + out = nifti_strdup("''") ; return out ; /* string?? */ } lout = 4 ; /* initialize length of output */ @@ -1837,6 +5297,10 @@ } } out = (char *)calloc(1,lout) ; /* allocate output string */ + if( !out ){ + fprintf(stderr,"** escapize_string: failed to alloc %d bytes\n",lout); + return NULL; + } out[0] = '\'' ; /* opening quote mark */ for( ii=0,jj=1 ; ii < lstr ; ii++ ){ switch( str[ii] ){ @@ -1861,26 +5325,30 @@ } /*---------------------------------------------------------------------------*/ -/* Dump the information in a NIFTI image header to an XML-ish ASCII string +/*! Dump the information in a NIFTI image header to an XML-ish ASCII string that can later be converted back into a NIFTI header in - nifti_image_from_ascii(). The resulting string can be free()-ed when - you are done with it. ------------------------------------------------------------------------------*/ - -char *nifti_image_to_ascii( nifti_image *nim ) + nifti_image_from_ascii(). + + The resulting string can be free()-ed when you are done with it. +*//*-------------------------------------------------------------------------*/ +char *nifti_image_to_ascii( const nifti_image *nim ) { char *buf , *ebuf ; int nbuf ; if( nim == NULL ) return NULL ; /* stupid caller */ - buf = (char *)calloc(1,65534); nbuf = 0; /* longer than needed, to be safe */ + buf = (char *)calloc(1,65534); nbuf = 0; /* longer than needed, to be safe */ + if( !buf ){ + fprintf(stderr,"** NITA: failed to alloc %d bytes\n",65534); + return NULL; + } sprintf( buf , "<nifti_image\n" ) ; /* XML-ish opener */ sprintf( buf+strlen(buf) , " nifti_type = '%s'\n" , - (nim->nifti_type == 1) ? "NIFTI-1+" - :(nim->nifti_type == 2) ? "NIFTI-1" - :(nim->nifti_type == 3) ? "NIFTI-1A" + (nim->nifti_type == NIFTI_FTYPE_NIFTI1_1) ? "NIFTI-1+" + :(nim->nifti_type == NIFTI_FTYPE_NIFTI1_2) ? "NIFTI-1" + :(nim->nifti_type == NIFTI_FTYPE_ASCII ) ? "NIFTI-1A" : "ANALYZE-7.5" ) ; /** Strings that we don't control (filenames, etc.) that might @@ -1888,7 +5356,7 @@ - A few special characters are replaced by XML-style escapes, using the function escapize_string(). - On input, function unescape_string() reverses this process. - - The result is that the NIFTI ASCII-format header is XML-compliant. **/ + - The result is that the NIFTI ASCII-format header is XML-compliant. */ ebuf = escapize_string(nim->fname) ; sprintf( buf+strlen(buf) , " header_filename = %s\n",ebuf); free(ebuf); @@ -2031,10 +5499,10 @@ " qoffset_y = '%g'\n" " qoffset_z = '%g'\n" " qfac = '%g'\n" , - nim->quatern_b , nim->quatern_c , nim->quatern_c , + nim->quatern_b , nim->quatern_c , nim->quatern_d , nim->qoffset_x , nim->qoffset_y , nim->qoffset_z , nim->qfac ) ; - mat44_to_orientation( nim->qto_xyz , &i,&j,&k ) ; + nifti_mat44_to_orientation( nim->qto_xyz , &i,&j,&k ) ; if( i > 0 && j > 0 && k > 0 ) sprintf( buf+strlen(buf) , " qform_i_orientation = '%s'\n" @@ -2073,7 +5541,7 @@ nim->sto_ijk.m[3][0] , nim->sto_ijk.m[3][1] , nim->sto_ijk.m[3][2] , nim->sto_ijk.m[3][3] ) ; - mat44_to_orientation( nim->sto_xyz , &i,&j,&k ) ; + nifti_mat44_to_orientation( nim->sto_xyz , &i,&j,&k ) ; if( i > 0 && j > 0 && k > 0 ) sprintf( buf+strlen(buf) , " sform_i_orientation = '%s'\n" @@ -2084,16 +5552,25 @@ nifti_orientation_string(k) ) ; } + sprintf( buf+strlen(buf) , " num_ext = '%d'\n", nim->num_ext ) ; + sprintf( buf+strlen(buf) , "/>\n" ) ; /* XML-ish closer */ nbuf = strlen(buf) ; buf = (char *)realloc((void *)buf, nbuf+1); /* cut back to proper length */ + if( !buf ) fprintf(stderr,"** NITA: failed to realloc %d bytes\n",nbuf+1); return buf ; } /*---------------------------------------------------------------------------*/ -int short_order(void) /* determine this CPU's byte order */ +/*----------------------------------------------------------------------*/ +/*! get the byte order for this CPU + + - LSB_FIRST means least significant byte, first (little endian) + - MSB_FIRST means most significant byte, first (big endian) +*//*--------------------------------------------------------------------*/ +int nifti_short_order(void) /* determine this CPU's byte order */ { union { unsigned char bb[2] ; short ss ; } fred ; @@ -2122,18 +5599,18 @@ put rhs string into nim->"nam" string, with max length = "ml" */ #define QSTR(nam,ml) if( strcmp(lhs,#nam) == 0 ) \ - strncpy(nim->nam,rhs,ml), nim->intent_name[ml]='\0' + strncpy(nim->nam,rhs,ml), nim->nam[ml]='\0' /*---------------------------------------------------------------------------*/ -/* Take an XML-ish ASCII string and create a NIFTI image header to match. - NULL is returned if enough information isn't present in the input string. +/*! Take an XML-ish ASCII string and create a NIFTI image header to match. + + NULL is returned if enough information isn't present in the input string. - The image data can later be loaded with nifti_image_load(). - The struct returned here can be liberated with nifti_image_free(). - Not a lot of error checking is done here to make sure that the input values are reasonable! ------------------------------------------------------------------------------*/ - -nifti_image *nifti_image_from_ascii( char *str ) +*//*-------------------------------------------------------------------------*/ +nifti_image *nifti_image_from_ascii( const char *str, int * bytes_read ) { char lhs[1024] , rhs[1024] ; int ii , spos, nn , slen ; @@ -2150,13 +5627,17 @@ /* create empty image struct */ nim = (nifti_image *) calloc( 1 , sizeof(nifti_image) ) ; + if( !nim ){ + fprintf(stderr,"** NIFA: failed to alloc nifti_image\n"); + return NULL; + } nim->nx = nim->ny = nim->nz = nim->nt = nim->nu = nim->nv = nim->nw = 1 ; nim->dx = nim->dy = nim->dz = nim->dt = nim->du = nim->dv = nim->dw = nim->qfac = 1.0 ; - nim->byteorder = short_order() ; + nim->byteorder = nifti_short_order() ; /* starting at str[spos], scan for "equations" of the form lhs = 'rhs' @@ -2196,16 +5677,20 @@ Start with special cases that don't fit the QNUM/QSTR macros. */ if( strcmp(lhs,"nifti_type") == 0 ){ - if( strcmp(rhs,"ANALYZE-7.5") == 0 ) nim->nifti_type = 0 ; - else if( strcmp(rhs,"NIFTI-1+") == 0 ) nim->nifti_type = 1 ; - else if( strcmp(rhs,"NIFTI-1") == 0 ) nim->nifti_type = 2 ; - else if( strcmp(rhs,"NIFTI-1A") == 0 ) nim->nifti_type = 3 ; + if( strcmp(rhs,"ANALYZE-7.5") == 0 ) + nim->nifti_type = NIFTI_FTYPE_ANALYZE ; + else if( strcmp(rhs,"NIFTI-1+") == 0 ) + nim->nifti_type = NIFTI_FTYPE_NIFTI1_1 ; + else if( strcmp(rhs,"NIFTI-1") == 0 ) + nim->nifti_type = NIFTI_FTYPE_NIFTI1_2 ; + else if( strcmp(rhs,"NIFTI-1A") == 0 ) + nim->nifti_type = NIFTI_FTYPE_ASCII ; } else if( strcmp(lhs,"header_filename") == 0 ){ - nim->fname = strdup(rhs) ; + nim->fname = nifti_strdup(rhs) ; } else if( strcmp(lhs,"image_filename") == 0 ){ - nim->iname = strdup(rhs) ; + nim->iname = nifti_strdup(rhs) ; } else if( strcmp(lhs,"sto_xyz_matrix") == 0 ){ sscanf( rhs , "%f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f" , @@ -2269,15 +5754,18 @@ else QNUM(slice_start) ; else QNUM(slice_end) ; else QNUM(slice_duration) ; + else QNUM(num_ext) ; } /* end of while loop */ - /** do miscellaneous checking and cleanup **/ - - if( nim->ndim <= 0 ){ nifti_image_free(nim); return NULL; } /** bad! **/ + if( bytes_read ) *bytes_read = spos+1; /* "process" last '\n' */ + + /* do miscellaneous checking and cleanup */ + + if( nim->ndim <= 0 ){ nifti_image_free(nim); return NULL; } /* bad! */ nifti_datatype_sizes( nim->datatype, &(nim->nbyper), &(nim->swapsize) ); - if( nim->nbyper == 0 ){ nifti_image_free(nim); return NULL; } /** bad! **/ + if( nim->nbyper == 0 ){ nifti_image_free(nim); return NULL; } /* bad! */ nim->dim[0] = nim->ndim ; nim->dim[1] = nim->nx ; nim->pixdim[1] = nim->dx ; @@ -2292,20 +5780,608 @@ * nim->nt * nim->nu * nim->nv * nim->nw ; if( nim->qform_code > 0 ) - nim->qto_xyz = quatern_to_mat44( - nim->quatern_b, nim->quatern_c, nim->quatern_c, + nim->qto_xyz = nifti_quatern_to_mat44( + nim->quatern_b, nim->quatern_c, nim->quatern_d, nim->qoffset_x, nim->qoffset_y, nim->qoffset_z, nim->dx , nim->dy , nim->dz , nim->qfac ) ; else - nim->qto_xyz = quatern_to_mat44( + nim->qto_xyz = nifti_quatern_to_mat44( 0.0 , 0.0 , 0.0 , 0.0 , 0.0 , 0.0 , nim->dx , nim->dy , nim->dz , 0.0 ) ; - nim->qto_ijk = mat44_inverse( nim->qto_xyz ) ; + + nim->qto_ijk = nifti_mat44_inverse( nim->qto_xyz ) ; if( nim->sform_code > 0 ) - nim->sto_ijk = mat44_inverse( nim->sto_xyz ) ; + nim->sto_ijk = nifti_mat44_inverse( nim->sto_xyz ) ; return nim ; } + + +/*---------------------------------------------------------------------------*/ +/*! validate the nifti_image + + \return 1 if the structure seems valid, otherwise 0 + + \sa nifti_nim_has_valid_dims, nifti_hdr_looks_good +*//*-------------------------------------------------------------------------*/ +int nifti_nim_is_valid(nifti_image * nim, int complain) +{ + int errs = 0; + + if( !nim ){ + fprintf(stderr,"** is_valid_nim: nim is NULL\n"); + return 0; + } + + if( g_opts.debug > 2 ) fprintf(stderr,"-d nim_is_valid check...\n"); + + /**- check that dim[] matches the individual values ndim, nx, ny, ... */ + if( ! nifti_nim_has_valid_dims(nim,complain) ){ + if( !complain ) return 0; + errs++; + } + + /* might check nbyper, pixdim, q/sforms, swapsize, nifti_type, ... */ + + /**- be explicit in return of 0 or 1 */ + if( errs > 0 ) return 0; + else return 1; +} + +/*---------------------------------------------------------------------------*/ +/*! validate nifti dimensions + + \return 1 if valid, 0 if not + + \sa nifti_nim_is_valid, nifti_hdr_looks_good + + rely on dim[] as the master +*//*-------------------------------------------------------------------------*/ +int nifti_nim_has_valid_dims(nifti_image * nim, int complain) +{ + int c, prod, errs = 0; + + /**- start with dim[0]: failure here is considered terminal */ + if( nim->dim[0] <= 0 || nim->dim[0] > 7 ){ + errs++; + if( complain ) + fprintf(stderr,"** NVd: dim[0] (%d) out of range [1,7]\n",nim->dim[0]); + return 0; + } + + /**- check whether ndim equals dim[0] */ + if( nim->ndim != nim->dim[0] ){ + errs++; + if( ! complain ) return 0; + fprintf(stderr,"** NVd: ndim != dim[0] (%d,%d)\n",nim->ndim,nim->dim[0]); + } + + /**- compare each dim[i] to the proper nx, ny, ... */ + if( ( (nim->dim[0] >= 1) && (nim->dim[1] != nim->nx) ) || + ( (nim->dim[0] >= 2) && (nim->dim[2] != nim->ny) ) || + ( (nim->dim[0] >= 3) && (nim->dim[3] != nim->nz) ) || + ( (nim->dim[0] >= 4) && (nim->dim[4] != nim->nt) ) || + ( (nim->dim[0] >= 5) && (nim->dim[5] != nim->nu) ) || + ( (nim->dim[0] >= 6) && (nim->dim[6] != nim->nv) ) || + ( (nim->dim[0] >= 7) && (nim->dim[7] != nim->nw) ) ){ + errs++; + if( !complain ) return 0; + fprintf(stderr,"** NVd mismatch: dims = %d,%d,%d,%d,%d,%d,%d\n" + " nxyz... = %d,%d,%d,%d,%d,%d,%d\n", + nim->dim[1], nim->dim[2], nim->dim[3], + nim->dim[4], nim->dim[5], nim->dim[6], nim->dim[7], + nim->nx, nim->ny, nim->nz, + nim->nt, nim->nu, nim->nv, nim->nw ); + } + + /**- check the dimensions, and that their product matches nvox */ + prod = 1; + for( c = 1; c <= nim->dim[0]; c++ ){ + if( nim->dim[c] > 0) + prod *= nim->dim[c]; + else if( nim->dim[c] <= 0 ){ + if( !complain ) return 0; + fprintf(stderr,"** NVd: dim[%d] (=%d) <= 0\n",c, nim->dim[c]); + errs++; + } + } + if( prod != nim->nvox ){ + if( ! complain ) return 0; + fprintf(stderr,"** NVd: nvox does not match dimension product (%d, %d)\n", + nim->nvox, prod); + errs++; + } + + /**- if debug, warn about any remaining dim that is neither 0, nor 1 */ + /* (values in dims above dim[0] are undefined, as reminded by Cinly + Ooi and Alle Meije Wink) 16 Nov 2005 [rickr] */ + if( g_opts.debug > 1 ) + for( c = nim->dim[0]+1; c <= 7; c++ ) + if( nim->dim[c] != 0 && nim->dim[c] != 1 ) + fprintf(stderr,"** NVd warning: dim[%d] = %d, but ndim = %d\n", + c, nim->dim[c], nim->dim[0]); + + if( g_opts.debug > 2 ) + fprintf(stderr,"-d nim_has_valid_dims check, errs = %d\n", errs); + + /**- return invalid or valid */ + if( errs > 0 ) return 0; + else return 1; +} + + +/*---------------------------------------------------------------------------*/ +/*! read a nifti image, collapsed across dimensions according to dims[8] <pre> + + This function may be used to read parts of a nifti dataset, such as + the time series for a single voxel, or perhaps a slice. It is similar + to nifti_image_load(), though the passed 'data' parameter is used for + returning the image, not nim->data. + + \param nim given nifti_image struct, corresponding to the data file + \param dims given list of dimensions (see below) + \param data pointer to data pointer (if *data is NULL, data will be + allocated, otherwise not) + + Here, dims is an array of 8 ints, similar to nim->dim[8]. While dims[0] + is unused at this point, the other indices specify which dimensions to + collapse (and at which index), and which not to collapse. If dims[i] is + set to -1, then that entire dimension will be read in, from index 0 to + index (nim->dim[i] - 1). If dims[i] >= 0, then only that index will be + read in (so dims[i] must also be < nim->dim[i]). + + Example: given nim->dim[8] = { 4, 64, 64, 21, 80, 1, 1, 1 } (4-D dataset) + + if dims[8] = { 0, 5, 4, 17, -1, -1, -1, -1 } + -> read time series for voxel i,j,k = 5,4,17 + + if dims[8] = { 0, -1, -1, -1, 17, -1, -1, -1 } + -> read single volume at time point 17 + + Example: given nim->dim[8] = { 6, 64, 64, 21, 80, 4, 3, 1 } (6-D dataset) + + if dims[8] = { 0, 5, 4, 17, -1, 2, 1, 0 } + -> read time series for the voxel i,j,k = 5,4,17, and dim 5,6 = 2,1 + + if dims[8] = { 0, 5, 4, -1, -1, 0, 0, 0 } + -> read time series for slice at i,j = 5,4, and dim 5,6,7 = 0,0,0 + (note that dims[7] is not relevant, but must be 0 or -1) + + If *data is NULL, then *data will be set as a pointer to new memory, + allocated here for the resulting collapsed image data. + + e.g. { int dims[8] = { 0, 5, 4, 17, -1, -1, -1, -1 }; + void * data = NULL; + ret_val = nifti_read_collapsed_image(nim, dims, &data); + if( ret_val > 0 ){ + process_time_series(data); + if( data != NULL ) free(data); + } + } + + NOTE: If *data is not NULL, then it will be assumed that it points to + valid memory, sufficient to hold the results. This is done for + speed and possibly repeated calls to this function. + + e.g. { int dims[8] = { 0, -1, -1, -1, -1, -1, -1, -1 }; + void * data = NULL; + for( zslice = 0; zslice < nzslices; zslice++ ){ + dims[3] = zslice; + ret_val = nifti_read_collapsed_image(nim, dims, &data); + if( ret_val > 0 ) process_slice(zslice, data); + } + if( data != NULL ) free(data); + } + + \return + - the total number of bytes read, or < 0 on failure + - the read and byte-swapped data, in 'data' </pre> + + \sa nifti_image_read, nifti_image_free, nifti_image_read_bricks + nifti_image_load +*//*-------------------------------------------------------------------------*/ +int nifti_read_collapsed_image( nifti_image * nim, const int dims [8], + void ** data ) +{ + znzFile fp; + int pivots[8], prods[8], nprods; /* sizes are bounded by dims[], so 8 */ + int c, bytes; + + /** - check pointers for sanity */ + if( !nim || !dims || !data ){ + fprintf(stderr,"** nifti_RCI: bad params %p, %p, %p\n", + (void *)nim, (void *)dims, (void *)data); + return -1; + } + + if( g_opts.debug > 2 ){ + fprintf(stderr,"-d read_collapsed_image:\n dims ="); + for(c = 0; c < 8; c++) fprintf(stderr," %3d", dims[c]); + fprintf(stderr,"\n nim->dims ="); + for(c = 0; c < 8; c++) fprintf(stderr," %3d", nim->dim[c]); + fputc('\n', stderr); + } + + /** - verify that dim[] makes sense */ + if( ! nifti_nim_is_valid(nim, g_opts.debug > 0) ){ + fprintf(stderr,"** invalid nim (file is '%s')\n", nim->fname ); + return -1; + } + + /** - verify that dims[] makes sense for this dataset */ + for( c = 1; c <= nim->dim[0]; c++ ){ + if( dims[c] >= nim->dim[c] ){ + fprintf(stderr,"** nifti_RCI: dims[%d] >= nim->dim[%d] (%d,%d)\n", + c, c, dims[c], nim->dim[c]); + return -1; + } + } + + /** - prepare pivot list - pivots are fixed indices */ + if( make_pivot_list(nim, dims, pivots, prods, &nprods) < 0 ) return -1; + + bytes = rci_alloc_mem(data, prods, nprods, nim->nbyper); + if( bytes < 0 ) return -1; + + /** - open the image file for reading at the appropriate offset */ + fp = nifti_image_load_prep( nim ); + if( ! fp ){ free(*data); *data = NULL; return -1; } /* failure */ + + /** - call the recursive reading function, passing nim, the pivot info, + location to store memory, and file pointer and position */ + c = rci_read_data(nim, pivots,prods,nprods,dims, + (char *)*data, fp, znztell(fp)); + znzclose(fp); /* in any case, close the file */ + if( c < 0 ){ free(*data); *data = NULL; return -1; } /* failure */ + + if( g_opts.debug > 1 ) + fprintf(stderr,"+d read %d bytes of collapsed image from %s\n", + bytes, nim->fname); + + return bytes; +} + + +/* read the data from the file pointed to by fp + + - this a recursive function, so start with the base case + - data is now (char *) for easy incrementing + + return 0 on success, < 0 on failure +*/ +static int rci_read_data(nifti_image * nim, int * pivots, int * prods, + int nprods, const int dims[], char * data, znzFile fp, int base_offset) +{ + int c, sublen, offset, read_size; + + /* bad check first - base_offset may not have been checked */ + if( base_offset < 0 || nprods <= 0 ){ + fprintf(stderr,"** rci_read_data, bad params, %d,%d\n", + nprods, base_offset); + return -1; + } + + /* base case: actually read the data */ + if( nprods == 1 ){ + int nread, bytes; + + /* make sure things look good here */ + if( *pivots != 0 ){ + fprintf(stderr,"** rciRD: final pivot == %d!\n", *pivots); + return -1; + } + + /* so just seek and read (prods[0] * nbyper) bytes from the file */ + znzseek(fp, base_offset, SEEK_SET); + bytes = prods[0] * nim->nbyper; + nread = nifti_read_buffer(fp, data, bytes, nim); + if( nread != bytes ){ + fprintf(stderr,"** rciRD: read only %d of %d bytes from '%s'\n", + nread, bytes, nim->fname); + return -1; + } else if( g_opts.debug > 3 ) + fprintf(stderr,"+d successful read of %d bytes at offset %d\n", + bytes, base_offset); + + return 0; /* done with base case - return success */ + } + + /* not the base case, so do a set of reduced reads */ + + /* compute size of sub-brick: all dimensions below pivot */ + for( c = 1, sublen = 1; c < *pivots; c++ ) sublen *= nim->dim[c]; + + /* compute number of values to read, i.e. remaining prods */ + for( c = 1, read_size = 1; c < nprods; c++ ) read_size *= prods[c]; + read_size *= nim->nbyper; /* and multiply by bytes per voxel */ + + /* now repeatedly compute offsets, and recursively read */ + for( c = 0; c < prods[0]; c++ ){ + /* offset is (c * sub-block size (including pivot dim)) */ + /* + (dims[] index into pivot sub-block) */ + /* the unneeded multiplication is to make this more clear */ + offset = c * sublen * nim->dim[*pivots] + sublen * dims[*pivots]; + offset *= nim->nbyper; + + if( g_opts.debug > 3 ) + fprintf(stderr,"-d reading %d bytes, foff %d + %d, doff %d\n", + read_size, base_offset, offset, c*read_size); + + /* now read the next level down, adding this offset */ + if( rci_read_data(nim, pivots+1, prods+1, nprods-1, dims, + data + c * read_size, fp, base_offset + offset) < 0 ) + return -1; + } + + return 0; +} + + +/* allocate memory for all collapsed image data + + If *data is already set, do not allocate, but still calculate + size for debug report. + + return total size on success, and < 0 on failure +*/ +static int rci_alloc_mem(void ** data, int prods[8], int nprods, int nbyper ) +{ + int size, index; + + if( nbyper < 0 || nprods < 1 || nprods > 8 ){ + fprintf(stderr,"** rci_am: bad params, %d, %d\n", nbyper, nprods); + return -1; + } + + for( index = 0, size = 1; index < nprods; index++ ) + size *= prods[index]; + + size *= nbyper; + + if( ! *data ){ /* then allocate what is needed */ + if( g_opts.debug > 1 ) + fprintf(stderr,"+d alloc %d (= %d x %d) bytes for collapsed image\n", + size, size/nbyper, nbyper); + + *data = malloc(size); /* actually allocate the memory */ + if( ! *data ){ + fprintf(stderr,"** rci_am: failed to alloc %d bytes for data\n", size); + return -1; + } + } else if( g_opts.debug > 1 ) + fprintf(stderr,"-d rci_am: *data already set, need %d (%d x %d) bytes\n", + size, size/nbyper, nbyper); + + return size; +} + + +/* prepare a pivot list for reading + + The pivot points are the indices into dims where the calling function + wants to collapse a dimension. The last pivot should always be zero + (note that we have space for that in the lists). +*/ +static int make_pivot_list(nifti_image * nim, const int dims[], int pivots[], + int prods[], int * nprods ) +{ + int len, index; + + len = 0; + index = nim->dim[0]; + while( index > 0 ){ + prods[len] = 1; + while( index > 0 && (nim->dim[index] == 1 || dims[index] == -1) ){ + prods[len] *= nim->dim[index]; + index--; + } + pivots[len] = index; + len++; + index--; /* fine, let it drop out at -1 */ + } + + /* make sure to include 0 as a pivot (instead of just 1, if it is) */ + if( pivots[len-1] != 0 ){ + pivots[len] = 0; + prods[len] = 1; + len++; + } + + *nprods = len; + + if( g_opts.debug > 2 ){ + fprintf(stderr,"+d pivot list created, pivots :"); + for(index = 0; index < len; index++) fprintf(stderr," %d", pivots[index]); + fprintf(stderr,", prods :"); + for(index = 0; index < len; index++) fprintf(stderr," %d", prods[index]); + fputc('\n',stderr); + } + + return 0; +} + + +#undef ISEND +#define ISEND(c) ( (c)==']' || (c)=='}' || (c)=='\0' ) + +/*---------------------------------------------------------------------*/ +/*! Get an integer list in the range 0..(nvals-1), from the + character string str. If we call the output pointer fred, + then fred[0] = number of integers in the list (> 0), and + fred[i] = i-th integer in the list for i=1..fred[0]. + If on return, fred == NULL or fred[0] == 0, then something is + wrong, and the caller must deal with that. + + Syntax of input string: + - initial '{' or '[' is skipped, if present + - ends when '}' or ']' or end of string is found + - contains entries separated by commas + - entries have one of these forms: + - a single number + - a dollar sign '$', which means nvals-1 + - a sequence of consecutive numbers in the form "a..b" or + "a-b", where "a" and "b" are single numbers (or '$') + - a sequence of evenly spaced numbers in the form + "a..b(c)" or "a-b(c)", where "c" encodes the step + - Example: "[2,7..4,3..9(2)]" decodes to the list + 2 7 6 5 4 3 5 7 9 + - entries should be in the range 0..nvals-1 + + (borrowed, with permission, from thd_intlist.c) +*//*-------------------------------------------------------------------*/ +int * nifti_get_intlist( int nvals , const char * str ) +{ + int *subv = NULL ; + int ii , ipos , nout , slen ; + int ibot,itop,istep , nused ; + char *cpt ; + + /* Meaningless input? */ + if( nvals < 1 ) return NULL ; + + /* No selection list? */ + if( str == NULL || str[0] == '\0' ) return NULL ; + + /* skip initial '[' or '{' */ + subv = (int *) malloc( sizeof(int) * 2 ) ; + subv[0] = nout = 0 ; + + ipos = 0 ; + if( str[ipos] == '[' || str[ipos] == '{' ) ipos++ ; + + if( g_opts.debug > 1 ) + fprintf(stderr,"-d making int_list (vals = %d) from '%s'\n", nvals, str); + + /**- for each sub-selector until end of input... */ + + slen = strlen(str) ; + while( ipos < slen && !ISEND(str[ipos]) ){ + + while( isspace(str[ipos]) ) ipos++ ; /* skip blanks */ + if( ISEND(str[ipos]) ) break ; /* done */ + + /**- get starting value */ + + if( str[ipos] == '$' ){ /* special case */ + ibot = nvals-1 ; ipos++ ; + } else { /* decode an integer */ + ibot = strtol( str+ipos , &cpt , 10 ) ; + if( ibot < 0 ){ + fprintf(stderr,"** ERROR: list index %d is out of range 0..%d\n", + ibot,nvals-1) ; + free(subv) ; return NULL ; + } + if( ibot >= nvals ){ + fprintf(stderr,"** ERROR: list index %d is out of range 0..%d\n", + ibot,nvals-1) ; + free(subv) ; return NULL ; + } + nused = (cpt-(str+ipos)) ; + if( ibot == 0 && nused == 0 ){ + fprintf(stderr,"** ERROR: list syntax error '%s'\n",str+ipos) ; + free(subv) ; return NULL ; + } + ipos += nused ; + } + + while( isspace(str[ipos]) ) ipos++ ; /* skip blanks */ + + /**- if that's it for this sub-selector, add one value to list */ + + if( str[ipos] == ',' || ISEND(str[ipos]) ){ + nout++ ; + subv = (int *) realloc( (char *)subv , sizeof(int) * (nout+1) ) ; + subv[0] = nout ; + subv[nout] = ibot ; + if( ISEND(str[ipos]) ) break ; /* done */ + ipos++ ; continue ; /* re-start loop at next sub-selector */ + } + + /**- otherwise, must have '..' or '-' as next inputs */ + + if( str[ipos] == '-' ){ + ipos++ ; + } else if( str[ipos] == '.' && str[ipos+1] == '.' ){ + ipos++ ; ipos++ ; + } else { + fprintf(stderr,"** ERROR: index list syntax is bad: '%s'\n", + str+ipos) ; + free(subv) ; return NULL ; + } + + /**- get ending value for loop now */ + + if( str[ipos] == '$' ){ /* special case */ + itop = nvals-1 ; ipos++ ; + } else { /* decode an integer */ + itop = strtol( str+ipos , &cpt , 10 ) ; + if( itop < 0 ){ + fprintf(stderr,"** ERROR: index %d is out of range 0..%d\n", + itop,nvals-1) ; + free(subv) ; return NULL ; + } + if( itop >= nvals ){ + fprintf(stderr,"** ERROR: index %d is out of range 0..%d\n", + itop,nvals-1) ; + free(subv) ; return NULL ; + } + nused = (cpt-(str+ipos)) ; + if( itop == 0 && nused == 0 ){ + fprintf(stderr,"** ERROR: index list syntax error '%s'\n",str+ipos) ; + free(subv) ; return NULL ; + } + ipos += nused ; + } + + /**- set default loop step */ + + istep = (ibot <= itop) ? 1 : -1 ; + + while( isspace(str[ipos]) ) ipos++ ; /* skip blanks */ + + /**- check if we have a non-default loop step */ + + if( str[ipos] == '(' ){ /* decode an integer */ + ipos++ ; + istep = strtol( str+ipos , &cpt , 10 ) ; + if( istep == 0 ){ + fprintf(stderr,"** ERROR: index loop step is 0!\n") ; + free(subv) ; return NULL ; + } + nused = (cpt-(str+ipos)) ; + ipos += nused ; + if( str[ipos] == ')' ) ipos++ ; + if( (ibot-itop)*istep > 0 ){ + fprintf(stderr,"** WARNING: index list '%d..%d(%d)' means nothing\n", + ibot,itop,istep ) ; + } + } + + /**- add values to output */ + + for( ii=ibot ; (ii-itop)*istep <= 0 ; ii += istep ){ + nout++ ; + subv = (int *) realloc( (char *)subv , sizeof(int) * (nout+1) ) ; + subv[0] = nout ; + subv[nout] = ii ; + } + + /**- check if we have a comma to skip over */ + + while( isspace(str[ipos]) ) ipos++ ; /* skip blanks */ + if( str[ipos] == ',' ) ipos++ ; /* skip commas */ + + } /* end of loop through selector string */ + + if( g_opts.debug > 1 ) { + fprintf(stderr,"+d int_list (vals = %d): ", subv[0]); + for( ii = 1; ii <= subv[0]; ii++ ) fprintf(stderr,"%d ", subv[ii]); + fputc('\n',stderr); + } + + if( subv[0] == 0 ){ free(subv); subv = NULL; } + return subv ; +}
--- a/conversion/nifti1/nifti1_io.h +++ b/conversion/nifti1/nifti1_io.h @@ -1,3 +1,8 @@ +/** \file nifti1_io.h + \brief Data structures for using nifti1_io API. + - Written by Bob Cox, SSCC NIMH + - Revisions by Rick Reynolds, SSCC NIMH + */ #ifndef _NIFTI_IO_HEADER_ #define _NIFTI_IO_HEADER_ @@ -12,6 +17,8 @@ #endif #include "nifti1.h" /*** NIFTI-1 header specification ***/ +#include "znzlib.h" + /*=================*/ #ifdef __cplusplus extern "C" { @@ -32,6 +39,20 @@ /***** incidental or otherwise, caused by any use of this document. *****/ /*****===================================================================*****/ +/* + Modified by: Mark Jenkinson (FMRIB Centre, University of Oxford, UK) + Date: July/August 2004 + + Mainly adding low-level IO and changing things to allow gzipped files + to be read and written + Full backwards compatability should have been maintained + + Modified by: Rick Reynolds (SSCC/DIRP/NIMH, National Institutes of Health) + Date: December 2004 + + Modified and added many routines for I/O. +*/ + /********************** Some sample data structures **************************/ typedef struct { /** 4x4 matrix struct **/ @@ -44,65 +65,105 @@ /*...........................................................................*/ -typedef struct { /** Image storage struct **/ +/*! \struct nifti_image + \brief High level data structure for open nifti datasets in the + nifti1_io API. Note that this structure is not part of the + nifti1 format definition; it is used to implement one API + for reading/writing formats in the nifti1 format. + */ +typedef struct { /*!< Image storage struct **/ - int ndim ; /* last dimension greater than 1 (1..7) */ - int nx,ny,nz , nt,nu,nv,nw ; /* dimensions of grid array */ - int dim[8] ; /* dim[0]=ndim, dim[1]=nx, etc. */ - int nvox ; /* number of voxels = nx*ny*nz*... */ - int nbyper ; /* bytes per voxel */ - int datatype ; /* type of data in voxels: DT_* code */ - - float dx,dy,dz , dt,du,dv,dw ; /* grid spacings */ - float pixdim[8] ; /* pixdim[1]=dx, etc. */ + int ndim ; /*!< last dimension greater than 1 (1..7) */ + int nx ; /*!< dimensions of grid array */ + int ny ; /*!< dimensions of grid array */ + int nz ; /*!< dimensions of grid array */ + int nt ; /*!< dimensions of grid array */ + int nu ; /*!< dimensions of grid array */ + int nv ; /*!< dimensions of grid array */ + int nw ; /*!< dimensions of grid array */ + int dim[8] ; /*!< dim[0]=ndim, dim[1]=nx, etc. */ + int nvox ; /*!< number of voxels = nx*ny*nz*...*nw */ + int nbyper ; /*!< bytes per voxel, matches datatype */ + int datatype ; /*!< type of data in voxels: DT_* code */ - float scl_slope , scl_inter ; /* scaling parameters */ - - float cal_min , cal_max ; /* calibration parameters */ - - int qform_code , sform_code ; /* codes for (x,y,z) space meaning */ + float dx ; /*!< grid spacings */ + float dy ; /*!< grid spacings */ + float dz ; /*!< grid spacings */ + float dt ; /*!< grid spacings */ + float du ; /*!< grid spacings */ + float dv ; /*!< grid spacings */ + float dw ; /*!< grid spacings */ + float pixdim[8] ; /*!< pixdim[1]=dx, etc. */ - int freq_dim , /* indexes (1,2,3, or 0) for MRI */ - phase_dim , slice_dim ; /* directions in dim[]/pixdim[] */ + float scl_slope ; /*!< scaling parameter - slope */ + float scl_inter ; /*!< scaling parameter - intercept */ + + float cal_min ; /*!< calibration parameter, minimum */ + float cal_max ; /*!< calibration parameter, maximum */ - int slice_code ; /* code for slice timing pattern */ - int slice_start , slice_end ; /* indexes for start & stop of slices */ - float slice_duration ; /* time between individual slices */ + int qform_code ; /*!< codes for (x,y,z) space meaning */ + int sform_code ; /*!< codes for (x,y,z) space meaning */ + + int freq_dim ; /*!< indexes (1,2,3, or 0) for MRI */ + int phase_dim ; /*!< directions in dim[]/pixdim[] */ + int slice_dim ; /*!< directions in dim[]/pixdim[] */ - float quatern_b, quatern_c, /* quaternion transform parameters */ - quatern_d, qoffset_x, /* [when writing a dataset, these ] */ - qoffset_y, qoffset_z, /* [are used for qform, NOT qto_xyz] */ - qfac ; + int slice_code ; /*!< code for slice timing pattern */ + int slice_start ; /*!< index for start of slices */ + int slice_end ; /*!< index for end of slices */ + float slice_duration ; /*!< time between individual slices */ + + /*! quaternion transform parameters + [when writing a dataset, these are used for qform, NOT qto_xyz] */ + float quatern_b , quatern_c , quatern_d , + qoffset_x , qoffset_y , qoffset_z , + qfac ; - mat44 qto_xyz ; /* qform: transform (i,j,k) to (x,y,z) */ - mat44 qto_ijk ; /* qform: transform (x,y,z) to (i,j,k) */ + mat44 qto_xyz ; /*!< qform: transform (i,j,k) to (x,y,z) */ + mat44 qto_ijk ; /*!< qform: transform (x,y,z) to (i,j,k) */ - mat44 sto_xyz ; /* sform: transform (i,j,k) to (x,y,z) */ - mat44 sto_ijk ; /* sform: transform (x,y,z) to (i,j,k) */ + mat44 sto_xyz ; /*!< sform: transform (i,j,k) to (x,y,z) */ + mat44 sto_ijk ; /*!< sform: transform (x,y,z) to (i,j,k) */ - float toffset ; /* time coordinate offset */ + float toffset ; /*!< time coordinate offset */ - int xyz_units , time_units ; /* dx,dy,dz & dt units: NIFTI_UNITS_* code */ + int xyz_units ; /*!< dx,dy,dz units: NIFTI_UNITS_* code */ + int time_units ; /*!< dt units: NIFTI_UNITS_* code */ - int nifti_type ; /* 0==ANALYZE, 2==NIFTI-1 (2 files), - 1==NIFTI-1 (1 file) , + int nifti_type ; /*!< 0==ANALYZE, 1==NIFTI-1 (1 file), + 2==NIFTI-1 (2 files), 3==NIFTI-ASCII (1 file) */ - int intent_code ; /* statistic type (or something) */ - float intent_p1, intent_p2, /* intent parameters */ - intent_p3 ; - char intent_name[16] ; + int intent_code ; /*!< statistic type (or something) */ + float intent_p1 ; /*!< intent parameters */ + float intent_p2 ; /*!< intent parameters */ + float intent_p3 ; /*!< intent parameters */ + char intent_name[16] ; /*!< optional description of intent data */ - char descrip[80], aux_file[24]; + char descrip[80] ; /*!< optional text to describe dataset */ + char aux_file[24] ; /*!< auxiliary filename */ - char *fname ; /* header filename (.hdr or .nii) */ - char *iname ; /* image filename (.img or .nii) */ - int iname_offset ; /* offset into iname where data starts */ - int swapsize ; /* swapping unit in image data (might be 0) */ - int byteorder ; /* byte order on disk (MSB_ or LSB_FIRST) */ - void *data ; /* pointer to data: nbyper*nvox bytes */ + char *fname ; /*!< header filename (.hdr or .nii) */ + char *iname ; /*!< image filename (.img or .nii) */ + int iname_offset ; /*!< offset into iname where data starts */ + int swapsize ; /*!< swap unit in image data (might be 0) */ + int byteorder ; /*!< byte order on disk (MSB_ or LSB_FIRST) */ + void *data ; /*!< pointer to data: nbyper*nvox bytes */ + + int num_ext ; /*!< number of extensions in ext_list */ + nifti1_extension * ext_list ; /*!< array of extension structs (with data) */ } nifti_image ; + + +/* struct for return from nifti_image_read_bricks() */ +typedef struct { + int nbricks; /* the number of allocated pointers in 'bricks' */ + int bsize; /* the length of each data block, in bytes */ + void ** bricks; /* array of pointers to data blocks */ +} nifti_brick_list; + + /*****************************************************************************/ /*--------------- Prototypes of functions defined in this file --------------*/ @@ -113,59 +174,111 @@ char *nifti_slice_string ( int ss ) ; char *nifti_orientation_string( int ii ) ; -int nifti_is_inttype( int dt ) ; +int nifti_is_inttype( int dt ) ; -mat44 mat44_inverse( mat44 R ) ; +mat44 nifti_mat44_inverse( mat44 R ) ; -mat33 mat33_inverse( mat33 R ) ; -mat33 mat33_polar ( mat33 A ) ; -float mat33_rownorm( mat33 A ) ; -float mat33_colnorm( mat33 A ) ; -float mat33_determ ( mat33 R ) ; -mat33 mat33_mul ( mat33 A , mat33 B ) ; +mat33 nifti_mat33_inverse( mat33 R ) ; +mat33 nifti_mat33_polar ( mat33 A ) ; +float nifti_mat33_rownorm( mat33 A ) ; +float nifti_mat33_colnorm( mat33 A ) ; +float nifti_mat33_determ ( mat33 R ) ; +mat33 nifti_mat33_mul ( mat33 A , mat33 B ) ; -void swap_2bytes ( int n , void *ar ) ; -void swap_4bytes ( int n , void *ar ) ; -void swap_8bytes ( int n , void *ar ) ; -void swap_16bytes( int n , void *ar ) ; -void swap_Nbytes ( int n , int siz , void *ar ) ; +void nifti_swap_2bytes ( int n , void *ar ) ; +void nifti_swap_4bytes ( int n , void *ar ) ; +void nifti_swap_8bytes ( int n , void *ar ) ; +void nifti_swap_16bytes( int n , void *ar ) ; +void nifti_swap_Nbytes ( int n , int siz , void *ar ) ; + +void swap_nifti_header ( struct nifti_1_header *h , int is_nifti ) ; +int nifti_get_filesize( const char *pathname ) ; -void swap_nifti_header( struct nifti_1_header *h , int is_nifti ) ; -unsigned int get_filesize( char *pathname ) ; +/* main read/write routines */ -nifti_image *nifti_image_read ( char *hname , int read_data ) ; -void nifti_image_load ( nifti_image *nim ) ; +nifti_image *nifti_image_read_bricks(const char *hname , int nbricks, + const int *blist, nifti_brick_list * NBL); +int nifti_image_load_bricks(nifti_image *nim , int nbricks, + const int *blist, nifti_brick_list * NBL); +void nifti_free_NBL( nifti_brick_list * NBL ); + +nifti_image *nifti_image_read ( const char *hname , int read_data ) ; +int nifti_image_load ( nifti_image *nim ) ; void nifti_image_unload ( nifti_image *nim ) ; void nifti_image_free ( nifti_image *nim ) ; -void nifti_image_write ( nifti_image *nim ) ; -void nifti_image_infodump( nifti_image *nim ) ; + +int nifti_read_collapsed_image( nifti_image * nim, const int dims [8], + void ** data ); + +void nifti_image_write ( nifti_image * nim ) ; +void nifti_image_write_bricks(nifti_image * nim, + const nifti_brick_list * NBL); +void nifti_image_infodump( const nifti_image * nim ) ; + +void nifti_disp_lib_hist( void ) ; /* to display library history */ +void nifti_disp_lib_version( void ) ; /* to display library version */ +int nifti_disp_matrix_orient( const char * mesg, mat44 mat ); + +char * nifti_image_to_ascii ( const nifti_image * nim ) ; +nifti_image *nifti_image_from_ascii( const char * str, int * bytes_read ) ; + +size_t nifti_get_volsize(const nifti_image *nim) ; -char * nifti_image_to_ascii ( nifti_image *nim ) ; -nifti_image *nifti_image_from_ascii( char *str ) ; +/* basic file operations */ +int nifti_set_filenames(nifti_image * nim, const char * prefix, int check, + int set_byte_order); +char * nifti_makehdrname (const char * prefix, int nifti_type, int check, + int comp); +char * nifti_makeimgname (const char * prefix, int nifti_type, int check, + int comp); +int is_nifti_file (const char *hname); +char * nifti_find_file_extension(const char * name); +int nifti_is_complete_filename(const char* fname); +int nifti_validfilename(const char* fname); + +int disp_nifti_1_header(const char * info, const nifti_1_header * hp ) ; +void nifti_set_debug_level( int level ) ; +void nifti_set_skip_blank_ext( int skip ) ; -int is_nifti_file( char *hname ) ; +int valid_nifti_brick_list(nifti_image * nim , int nbricks, + const int * blist, int disp_error); + +/* znzFile operations */ +znzFile nifti_image_open(const char * hname, char * opts, nifti_image ** nim); +znzFile nifti_image_write_hdr_img(nifti_image *nim, int write_data, + const char* opts); +znzFile nifti_image_write_hdr_img2( nifti_image *nim , int write_opts , + const char* opts, znzFile imgfile, const nifti_brick_list * NBL); +size_t nifti_read_buffer(znzFile fp, void* datatptr, size_t ntot, + nifti_image *nim); +int nifti_write_all_data(znzFile fp, nifti_image * nim, + const nifti_brick_list * NBL); +size_t nifti_write_buffer(znzFile fp, const void * buffer, size_t numbytes); +nifti_image *nifti_read_ascii_image(znzFile fp, char *fname, int flen, + int read_data); +znzFile nifti_write_ascii_image(nifti_image *nim, const nifti_brick_list * NBL, + const char * opts, int write_data, int leave_open); + void nifti_datatype_sizes( int datatype , int *nbyper, int *swapsize ) ; -void mat44_to_quatern( mat44 R , - float *qb, float *qc, float *qd, - float *qx, float *qy, float *qz, - float *dx, float *dy, float *dz, float *qfac ) ; +void nifti_mat44_to_quatern( mat44 R , + float *qb, float *qc, float *qd, + float *qx, float *qy, float *qz, + float *dx, float *dy, float *dz, float *qfac ) ; -mat44 quatern_to_mat44( float qb, float qc, float qd, - float qx, float qy, float qz, - float dx, float dy, float dz, float qfac ); +mat44 nifti_quatern_to_mat44( float qb, float qc, float qd, + float qx, float qy, float qz, + float dx, float dy, float dz, float qfac ); -mat44 make_orthog_mat44( float r11, float r12, float r13 , - float r21, float r22, float r23 , - float r31, float r32, float r33 ) ; +mat44 nifti_make_orthog_mat44( float r11, float r12, float r13 , + float r21, float r22, float r23 , + float r31, float r32, float r33 ) ; -int unescape_string( char *str ) ; /* string utility functions */ -char *escapize_string( char *str ) ; +int nifti_short_order(void) ; /* CPU byte order */ -int short_order(void) ; /* CPU byte order */ -/* Orientation codes that might be returned from mat44_to_orientation(). */ +/* Orientation codes that might be returned from nifti_mat44_to_orientation().*/ #define NIFTI_L2R 1 /* Left to Right */ #define NIFTI_R2L 2 /* Right to Left */ @@ -174,14 +287,98 @@ #define NIFTI_I2S 5 /* Inferior to Superior */ #define NIFTI_S2I 6 /* Superior to Inferior */ -void mat44_to_orientation( mat44 R , int *icod, int *jcod, int *kcod ) ; +void nifti_mat44_to_orientation( mat44 R , int *icod, int *jcod, int *kcod ) ; + +/*--------------------- Low level IO routines ------------------------------*/ + +char * nifti_findhdrname (const char* fname); +char * nifti_findimgname (const char* fname , int nifti_type); +int nifti_is_gzfile (const char* fname); + +char * nifti_makebasename(const char* fname); + + +/* other routines */ +struct nifti_1_header nifti_convert_nim2nhdr(const nifti_image* nim); +nifti_1_header * nifti_read_header(const char *hname, int *swapped, int check); +nifti_image * nifti_copy_nim_info(const nifti_image * src); +nifti_image * nifti_simple_init_nim(void); +nifti_image * nifti_convert_nhdr2nim(struct nifti_1_header nhdr, + const char * fname); + +int nifti_hdr_looks_good (const nifti_1_header * hdr); +int nifti_is_valid_ecode (int ecode); +int nifti_nim_is_valid (nifti_image * nim, int complain); +int nifti_nim_has_valid_dims (nifti_image * nim, int complain); +int is_valid_nifti_type (int nifti_type); +int nifti_type_and_names_match (nifti_image * nim, int show_warn); +int nifti_update_dims_from_array(nifti_image * nim); +void nifti_set_iname_offset (nifti_image *nim); +int nifti_set_type_from_names (nifti_image * nim); +int nifti_add_extension(nifti_image * nim, const char * data, int len, + int ecode ); +int nifti_copy_extensions (nifti_image *nim_dest,const nifti_image *nim_src); +int nifti_free_extensions (nifti_image *nim); +int * nifti_get_intlist (int nvals , const char *str); +char * nifti_strdup (const char *str); +int valid_nifti_extensions(const nifti_image *nim); + /*-------------------- Some C convenience macros ----------------------------*/ +/* NIfTI-1.1 extension codes: + see http://nifti.nimh.nih.gov/nifti-1/documentation/faq#Q21 */ + +#define NIFTI_ECODE_IGNORE 0 /* changed from UNKNOWN, 29 June 2005 */ + +#define NIFTI_ECODE_DICOM 2 /* intended for raw DICOM attributes */ + +#define NIFTI_ECODE_AFNI 4 /* Robert W Cox: rwcox@nih.gov + http://afni.nimh.nih.gov/afni */ + +#define NIFTI_ECODE_COMMENT 6 /* plain ASCII text only */ + +#define NIFTI_ECODE_XCEDE 8 /* David B Keator: dbkeator@uci.edu + http://www.nbirn.net/Resources + /Users/Applications/ + /xcede/index.htm */ + +#define NIFTI_ECODE_JIMDIMINFO 10 /* Mark A Horsfield: + mah5@leicester.ac.uk + http://someplace/something */ + +#define NIFTI_ECODE_WORKFLOW_FWDS 12 /* Kate Fissell: fissel+@pitt.edu + http://kraepelin.wpic.pitt.edu + /~fissell/NIFTI_ECODE_WORKFLOW_FWDS + /NIFTI_ECODE_WORKFLOW_FWDS.html */ + +#define NIFTI_MAX_ECODE 12 /******* maximum extension code *******/ + +/* nifti_type file codes */ +#define NIFTI_FTYPE_ANALYZE 0 +#define NIFTI_FTYPE_NIFTI1_1 1 +#define NIFTI_FTYPE_NIFTI1_2 2 +#define NIFTI_FTYPE_ASCII 3 +#define NIFTI_MAX_FTYPE 3 /* this should match the maximum code */ + +/*------------------------------------------------------------------------*/ +/*-- the rest of these apply only to nifti1_io.c, check for _NIFTI1_IO_C_ */ +/* Feb 9, 2005 [rickr] */ +#ifdef _NIFTI1_IO_C_ + +typedef struct { + int debug; /*!< debug level for status reports */ + int skip_blank_ext; /*!< skip extender if no extensions */ +} nifti_global_options; + +#undef LNI_FERR /* local nifti file error, to be compact and repetative */ +#define LNI_FERR(func,msg,file) \ + fprintf(stderr,"** ERROR (%s): %s '%s'\n",func,msg,file) + #undef swap_2 #undef swap_4 -#define swap_2(s) swap_2bytes(1,&(s)) /* s is a 2-byte short; swap in place */ -#define swap_4(v) swap_4bytes(1,&(v)) /* v is a 4-byte value; swap in place */ +#define swap_2(s) nifti_swap_2bytes(1,&(s)) /* s: 2-byte short; swap in place */ +#define swap_4(v) nifti_swap_4bytes(1,&(v)) /* v: 4-byte value; swap in place */ /***** isfinite() is a C99 macro, which is present in many C implementations already *****/ @@ -207,6 +404,11 @@ #define MSB_FIRST 2 #define REVERSE_ORDER(x) (3-(x)) /* convert MSB_FIRST <--> LSB_FIRST */ +#define LNI_MAX_NIA_EXT_LEN 100000 /* consider a longer extension invalid */ + +#endif /* _NIFTI1_IO_C_ section */ +/*------------------------------------------------------------------------*/ + /*=================*/ #ifdef __cplusplus }
--- a/conversion/nifti1/nii2mnc.c +++ b/conversion/nifti1/nii2mnc.c @@ -230,7 +230,7 @@ ss = ana_hdr.dime.dim[0]; if (ss != 0) { if (ss < 0 || ss > 7) { - swap_2(ss); + nifti_swap_2bytes(1, &(ss)); if (ss < 0 || ss > 7) { /* We should never get here!! */ fprintf(stderr, "Bad dimension count!!\n"); @@ -243,7 +243,7 @@ else { ss = ana_hdr.hk.sizeof_hdr; if (ss != sizeof(ana_hdr)) { - swap_4(ss); + nifti_swap_4bytes(1, &(ss)); if (ss != sizeof(ana_hdr)) { /* We should never get here!! */ fprintf(stderr, "Bad header size!!\n"); @@ -255,21 +255,21 @@ } if (must_swap) { - swap_4(ana_hdr.hk.sizeof_hdr) ; - swap_4(ana_hdr.hk.extents) ; - swap_2(ana_hdr.hk.session_error) ; + nifti_swap_4bytes(1, &(ana_hdr.hk.sizeof_hdr)); + nifti_swap_4bytes(1, &(ana_hdr.hk.extents)); + nifti_swap_2bytes(1, &(ana_hdr.hk.session_error)); - swap_4(ana_hdr.dime.compressed) ; - swap_4(ana_hdr.dime.verified) ; - swap_4(ana_hdr.dime.glmax) ; - swap_4(ana_hdr.dime.glmin) ; - swap_2bytes( 8 , ana_hdr.dime.dim ) ; - swap_4bytes( 8 , ana_hdr.dime.pixdim ) ; - swap_2(ana_hdr.dime.datatype) ; - swap_2(ana_hdr.dime.bitpix) ; - swap_4(ana_hdr.dime.vox_offset); - swap_4(ana_hdr.dime.cal_max); - swap_4(ana_hdr.dime.cal_min); + nifti_swap_4bytes(1, &(ana_hdr.dime.compressed)); + nifti_swap_4bytes(1, &(ana_hdr.dime.verified)); + nifti_swap_4bytes(1, &(ana_hdr.dime.glmax)); + nifti_swap_4bytes(1, &(ana_hdr.dime.glmin)); + nifti_swap_2bytes(8, ana_hdr.dime.dim); + nifti_swap_4bytes(8, ana_hdr.dime.pixdim); + nifti_swap_2bytes(1, &(ana_hdr.dime.datatype)); + nifti_swap_2bytes(1, &(ana_hdr.dime.bitpix)); + nifti_swap_4bytes(1, &(ana_hdr.dime.vox_offset)); + nifti_swap_4bytes(1, &(ana_hdr.dime.cal_max)); + nifti_swap_4bytes(1, &(ana_hdr.dime.cal_min)); } if (!qflag) {
new file mode 100644 --- /dev/null +++ b/conversion/nifti1/znzlib.c @@ -0,0 +1,274 @@ +/** \file znzlib.c + \brief Low level i/o interface to compressed and noncompressed files. + Written by Mark Jenkinson, FMRIB + +This library provides an interface to both compressed (gzip/zlib) and +uncompressed (normal) file IO. The functions are written to have the +same interface as the standard file IO functions. + +To use this library instead of normal file IO, the following changes +are required: + - replace all instances of FILE* with znzFile + - change the name of all function calls, replacing the initial character + f with the znz (e.g. fseek becomes znzseek) + one exception is rewind() -> znzrewind() + - add a third parameter to all calls to znzopen (previously fopen) + that specifies whether to use compression (1) or not (0) + - use znz_isnull rather than any (pointer == NULL) comparisons in the code + for znzfile types (normally done after a return from znzopen) + +NB: seeks for writable files with compression are quite restricted + + */ + +#include "znzlib.h" + +/* +znzlib.c (zipped or non-zipped library) + +***** This code is released to the public domain. ***** + +***** Author: Mark Jenkinson, FMRIB Centre, University of Oxford ***** +***** Date: September 2004 ***** + +***** Neither the FMRIB Centre, the University of Oxford, nor any of ***** +***** its employees imply any warranty of usefulness of this software ***** +***** for any purpose, and do not assume any liability for damages, ***** +***** incidental or otherwise, caused by any use of this document. ***** + +*/ + + +/* Note extra argument (use_compression) where + use_compression==0 is no compression + use_compression!=0 uses zlib (gzip) compression +*/ + +znzFile znzopen(const char *path, const char *mode, int use_compression) +{ + znzFile file; + file = (znzFile) calloc(1,sizeof(struct znzptr)); + if( file == NULL ){ + fprintf(stderr,"** ERROR: znzopen failed to alloc znzptr\n"); + return NULL; + } + + file->nzfptr = NULL; + +#ifdef HAVE_ZLIB + file->zfptr = NULL; + + if (use_compression) { + file->withz = 1; + if((file->zfptr = gzopen(path,mode)) == NULL) { + free(file); + file = NULL; + } + } else { +#endif + + file->withz = 0; + if((file->nzfptr = fopen(path,mode)) == NULL) { + free(file); + file = NULL; + } + +#ifdef HAVE_ZLIB + } +#endif + + return file; +} + + +znzFile znzdopen(int fd, const char *mode, int use_compression) +{ + znzFile file; + file = (znzFile) calloc(1,sizeof(struct znzptr)); + if( file == NULL ){ + fprintf(stderr,"** ERROR: znzdopen failed to alloc znzptr\n"); + return NULL; + } +#ifdef HAVE_ZLIB + if (use_compression) { + file->withz = 1; + file->zfptr = gzdopen(fd,mode); + file->nzfptr = NULL; + } else { +#endif + file->withz = 0; +#ifdef HAVE_FDOPEN + file->nzfptr = fdopen(fd,mode); +#endif +#ifdef HAVE_ZLIB + file->zfptr = NULL; + }; +#endif + return file; +} + + +int Xznzclose(znzFile * file) +{ + int retval = 0; + if (*file!=NULL) { +#ifdef HAVE_ZLIB + if ((*file)->zfptr!=NULL) { retval = gzclose((*file)->zfptr); } +#endif + if ((*file)->nzfptr!=NULL) { retval = fclose((*file)->nzfptr); } + + free(*file); + *file = NULL; + } + return retval; +} + + +size_t znzread(void* buf, size_t size, size_t nmemb, znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) + return (size_t) (gzread(file->zfptr,buf,((int) size)*((int) nmemb)) / size); +#endif + return fread(buf,size,nmemb,file->nzfptr); +} + +size_t znzwrite(const void* buf, size_t size, size_t nmemb, znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) + { + /* NOTE: We must typecast const away from the buffer because + gzwrite does not have complete const specification */ + return (size_t) ( gzwrite(file->zfptr,(void *)buf,size*nmemb) / size ); + } +#endif + return fwrite(buf,size,nmemb,file->nzfptr); +} + +long znzseek(znzFile file, long offset, int whence) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return (long) gzseek(file->zfptr,offset,whence); +#endif + return fseek(file->nzfptr,offset,whence); +} + +int znzrewind(znzFile stream) +{ + if (stream==NULL) { return 0; } +#ifdef HAVE_ZLIB + /* On some systems, gzrewind() fails for uncompressed files. + Use gzseek(), instead. 10, May 2005 [rickr] + + if (stream->zfptr!=NULL) return gzrewind(stream->zfptr); + */ + + if (stream->zfptr!=NULL) return (int)gzseek(stream->zfptr, 0L, SEEK_SET); +#endif + rewind(stream->nzfptr); + return 0; +} + +long znztell(znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return (long) gztell(file->zfptr); +#endif + return ftell(file->nzfptr); +} + +int znzputs(const char * str, znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzputs(file->zfptr,str); +#endif + return fputs(str,file->nzfptr); +} + + +char * znzgets(char* str, int size, znzFile file) +{ + if (file==NULL) { return NULL; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzgets(file->zfptr,str,size); +#endif + return fgets(str,size,file->nzfptr); +} + + +int znzflush(znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzflush(file->zfptr,Z_SYNC_FLUSH); +#endif + return fflush(file->nzfptr); +} + + +int znzeof(znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzeof(file->zfptr); +#endif + return feof(file->nzfptr); +} + + +int znzputc(int c, znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzputc(file->zfptr,c); +#endif + return fputc(c,file->nzfptr); +} + + +int znzgetc(znzFile file) +{ + if (file==NULL) { return 0; } +#ifdef HAVE_ZLIB + if (file->zfptr!=NULL) return gzgetc(file->zfptr); +#endif + return fgetc(file->nzfptr); +} + +#if !defined (WIN32) +int znzprintf(znzFile stream, const char *format, ...) +{ + int retval=0; + char *tmpstr; + va_list va; + if (stream==NULL) { return 0; } + va_start(va, format); +#ifdef HAVE_ZLIB + if (stream->zfptr!=NULL) { + int size; /* local to HAVE_ZLIB block */ + size = strlen(format) + 1000000; /* overkill I hope */ + tmpstr = (char *)calloc(1, size); + if( tmpstr == NULL ){ + fprintf(stderr,"** ERROR: znzprintf failed to alloc %d bytes\n", size); + return retval; + } + vsprintf(tmpstr,format,va); + retval=gzprintf(stream->zfptr,"%s",tmpstr); + free(tmpstr); + } else +#endif + { + retval=vfprintf(stream->nzfptr,format,va); + } + va_end(va); + return retval; +} + +#endif +
new file mode 100644 --- /dev/null +++ b/conversion/nifti1/znzlib.h @@ -0,0 +1,117 @@ +#ifndef _ZNZLIB_H_ +#define _ZNZLIB_H_ + +/* +znzlib.h (zipped or non-zipped library) + +***** This code is released to the public domain. ***** + +***** Author: Mark Jenkinson, FMRIB Centre, University of Oxford ***** +***** Date: September 2004 ***** + +***** Neither the FMRIB Centre, the University of Oxford, nor any of ***** +***** its employees imply any warranty of usefulness of this software ***** +***** for any purpose, and do not assume any liability for damages, ***** +***** incidental or otherwise, caused by any use of this document. ***** + +*/ + +/* + +This library provides an interface to both compressed (gzip/zlib) and +uncompressed (normal) file IO. The functions are written to have the +same interface as the standard file IO functions. + +To use this library instead of normal file IO, the following changes +are required: + - replace all instances of FILE* with znzFile + - change the name of all function calls, replacing the initial character + f with the znz (e.g. fseek becomes znzseek) + - add a third parameter to all calls to znzopen (previously fopen) + that specifies whether to use compression (1) or not (0) + - use znz_isnull rather than any (pointer == NULL) comparisons in the code + +NB: seeks for writable files with compression are quite restricted + +*/ + + +/*=================*/ +#ifdef __cplusplus +extern "C" { +#endif +/*=================*/ + +#include <stdio.h> +#include <stdlib.h> +#include <string.h> +#include <stdarg.h> + +#include "config.h" + +#ifdef HAVE_ZLIB +#if defined(ITKZLIB) +#include "itk_zlib.h" +#else +#include "zlib.h" +#endif +#endif + + +struct znzptr { + int withz; + FILE* nzfptr; +#ifdef HAVE_ZLIB + gzFile zfptr; +#endif +} ; + +/* the type for all file pointers */ +typedef struct znzptr * znzFile; + + +/* int znz_isnull(znzFile f); */ +/* int znzclose(znzFile f); */ +#define znz_isnull(f) ((f) == NULL) +#define znzclose(f) Xznzclose(&(f)) + +/* Note extra argument (use_compression) where + use_compression==0 is no compression + use_compression!=0 uses zlib (gzip) compression +*/ + +znzFile znzopen(const char *path, const char *mode, int use_compression); + +znzFile znzdopen(int fd, const char *mode, int use_compression); + +int Xznzclose(znzFile * file); + +size_t znzread(void* buf, size_t size, size_t nmemb, znzFile file); + +size_t znzwrite(const void* buf, size_t size, size_t nmemb, znzFile file); + +long znzseek(znzFile file, long offset, int whence); + +int znzrewind(znzFile stream); + +long znztell(znzFile file); + +int znzputs(const char *str, znzFile file); + +char * znzgets(char* str, int size, znzFile file); + +int znzputc(int c, znzFile file); + +int znzgetc(znzFile file); + +#if !defined(WIN32) +int znzprintf(znzFile stream, const char *format, ...); +#endif + +/*=================*/ +#ifdef __cplusplus +} +#endif +/*=================*/ + +#endif