Mercurial > hg > minc-tools
changeset 981:7c7801b48f7a
Initial revision
author | neelin <neelin> |
---|---|
date | Thu, 18 Jan 1996 14:52:14 +0000 |
parents | ba65debc47cf |
children | d993bc2a1e54 |
files | conversion/ecattominc/Makefile conversion/ecattominc/dump_ecat_header.c conversion/ecattominc/ecat_file.c conversion/ecattominc/ecat_file.h conversion/ecattominc/ecat_header_definition.h conversion/ecattominc/ecattominc.c conversion/ecattominc/insertblood.c |
diffstat | 7 files changed, 3294 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/Makefile @@ -0,0 +1,20 @@ +# -------------------------------------------------------------------- +# +# MINC Makefile +# + +ROOT = ../.. +include $(ROOT)/Make_machine_specific +include $(ROOT)/Make_configuration + +# Executable names +PROGS = ecattominc dump_ecat_header +EXTRA_OBJS = ecat_file.o insertblood.o +HEADERS = ecat_file.h ecat_header_definition.h +CDEFINES = -DDEBUG # cpp defines +LDOPT = $(PROG_LDOPT) +MANSECT = 1 +#MANPAGES = $(PROGS).$(MANSECT) + +include $(ROOT)/progs/Make_progs +
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/dump_ecat_header.c @@ -0,0 +1,150 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : dump_ecat_header +@INPUT : argc, argv - command line arguments +@OUTPUT : (none) +@RETURNS : error status +@DESCRIPTION: Dump ECAT file header information +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 10, 1996 (Peter Neelin) +@MODIFIED : $Log: dump_ecat_header.c,v $ +@MODIFIED : Revision 1.1 1996-01-18 14:52:14 neelin +@MODIFIED : Initial revision +@MODIFIED : +@COPYRIGHT : + Copyright 1996 Peter Neelin, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ + +#ifndef lint +static char rcsid[]="$Header: /private-cvsroot/minc/conversion/ecattominc/dump_ecat_header.c,v 1.1 1996-01-18 14:52:14 neelin Exp $"; +#endif + +#include <stdlib.h> +#include <stdio.h> +#include <string.h> +#include <errno.h> +#include <ecat_file.h> + +/* Main program */ + +int main(int argc, char *argv[]) +{ + Ecat_file *fp; + char *pname, *filename, *fieldname; + int index, volume, slice, iheader, ifield; + Ecat_field_name field; + char *description; + int length; + int error; + char svalue[ECAT_MAX_STRING_LENGTH]; + + /* Check the arguments */ + pname = argv[0]; + if ((argc < 2) || (argc > 6)) { + (void) fprintf(stderr, + "Usage: %s <ecat file> [<field> [<index> [<volume> <slice>]]]\n", + pname); + return EXIT_FAILURE; + } + filename = argv[1]; + fieldname = ((argc > 2) ? argv[2] : (char *) NULL); + index = ((argc > 3) ? atoi(argv[3]) : 0); + volume = ((argc > 4) ? atoi(argv[4]) : 0); + slice = ((argc > 5) ? atoi(argv[5]) : 0); + if ((fieldname != NULL) && (strlen(fieldname) == 0)) + fieldname = NULL; + + /* Open the file */ + fp = ecat_open(filename); + if (fp == NULL) { + (void) fprintf(stderr, "%s: Error opening file \"%s\"", + pname, filename); + if (errno != 0) { + (void) fprintf(stderr, " - "); + perror(NULL); + } + else { + (void) fprintf(stderr, "\n"); + } + exit(EXIT_FAILURE); + } + + /* Loop over main header then appropriate subheader */ + for (iheader = 0; iheader < 2; iheader++) { + + /* Print out header type */ + if (fieldname == NULL) { + if (iheader == 0) { + (void) printf("Main header:\n"); + } + else { + (void) printf("Subheader (volume %d, slice %d):\n", volume, slice); + } + } + + /* Loop over fields */ + ifield = -1; + do { + ifield++; + + /* Get information on the next field */ + if (iheader == 0) { + field = ecat_list_main(fp, ifield); + description = ecat_get_main_field_description(fp, field); + length = ecat_get_main_field_length(fp, field); + error = ecat_get_main_value(fp, field, index, + NULL, NULL, svalue); + } + else { + field = ecat_list_subhdr(fp, ifield); + description = ecat_get_subhdr_field_description(fp, field); + length = ecat_get_subhdr_field_length(fp, field); + error = ecat_get_subhdr_value(fp, volume, slice, field, index, + NULL, NULL, svalue); + } + + /* Check for the end of list */ + if (field == ECAT_No_Field) continue; + + /* If we are looking for a specific field, check to see if we've + found it */ + if ((fieldname != NULL) && (description != NULL) && + (strcmp(fieldname, description) != 0)) { + continue; + } + + /* Look for an error in reading the value */ + if (error) { + (void) fprintf(stderr, " Error retrieving field \"%s\".\n", + (description == NULL ? "<unknown>" : description)); + continue; + } + + /* Print out appropriate results */ + if ((fieldname == NULL) || (strcmp(fieldname, description) == 0)) { + (void) printf(" %s", description); + if (length > 1) + (void) printf("[%d/%d]", index, length); + (void) printf(" = %s\n", svalue); + } + + } while (field != ECAT_No_Field); + + } /* End of loop over headers */ + + /* Close the file */ + ecat_close(fp); + + return EXIT_SUCCESS; +} + +
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/ecat_file.c @@ -0,0 +1,1192 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_file.c +@DESCRIPTION: File containing routines to read ECAT image files +@GLOBALS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : $Log: ecat_file.c,v $ +@MODIFIED : Revision 1.1 1996-01-18 14:52:14 neelin +@MODIFIED : Initial revision +@MODIFIED : +@COPYRIGHT : + Copyright 1996 Peter Neelin, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ + +#ifndef lint +static char rcsid[]="$Header: /private-cvsroot/minc/conversion/ecattominc/ecat_file.c,v 1.1 1996-01-18 14:52:14 neelin Exp $"; +#endif + +#include <stdlib.h> +#include <stdio.h> +#include <string.h> +#include <memory.h> +#include <ecat_file.h> +#include <vax_conversions.h> + +/* Set some standard macros */ +#ifndef SEEK_SET +# define SEEK_SET 0 +#endif +#ifndef TRUE +# define TRUE 1 +# define FALSE 0 +#endif +#define MALLOC(size) ((void *) malloc_check(size)) +#define FREE(ptr) free(ptr) +#define REALLOC(ptr, size) ((void *) realloc_check(ptr, size)) +#define CALLOC(nelem, elsize) ((void *) calloc(nelem, elsize)) + +/* Constants */ +#define BLOCK_SIZE 512 +#define MAIN_HEADER_SIZE BLOCK_SIZE +#define SUBHEADER_SIZE BLOCK_SIZE +#define DIRBLOCK_SIZE BLOCK_SIZE +#define MAGIC_STRING "MATRIX" +#define FIRST_DIRBLOCK 2 +#define DRBLK_NEXT 1 +#define DRBLK_NUMUSED 3 +#define DRBLK_SKIP 4 +#define DRBLK_WIDTH 4 +#define DRBLK_PTR 1 + +/* Header definition types */ +typedef enum { + ecat_byte, ecat_short, ecat_long, ecat_float, ecat_char +} Ecat_type; + +typedef struct { + Ecat_field_name name; + int offset; + int length; + Ecat_type type; + char *description; +} Ecat_field_description_type; + +typedef struct { + int initialized; + int num_entries; + Ecat_field_description_type *fields; + Ecat_field_name *file_order; +} Ecat_header_table_type; + +typedef struct { + Ecat_header_table_type *main_header; + Ecat_header_table_type *subheader; +} Ecat_header_description_type; + +/* Header definition */ +#include "ecat_header_definition.h" + +/* ECAT file Type */ +struct Ecat_file { + FILE *file_pointer; + Ecat_header_description_type *header_description; + int vax_byte_order; + int num_planes; + int num_frames; + int num_bed_positions; + int num_gates; + int num_volumes; + unsigned char *main_header; + long cur_subhdr_offset; + unsigned char *subheader; + int num_subhdrs; + long *subhdr_offsets; +}; + +typedef enum { + ECAT_MAIN_HEADER, ECAT_SUBHEADER +} Ecat_which_header; + +/* Private functions */ +private Ecat_field_name ecat_list_fields(Ecat_file *file, + Ecat_which_header which_header, + int index); +private int ecat_get_field_length(Ecat_file *file, + Ecat_which_header which_header, + Ecat_field_name field); +private char *ecat_get_field_description(Ecat_file *file, + Ecat_which_header which_header, + Ecat_field_name field); +private int ecat_get_value(Ecat_file *file, + Ecat_which_header which_header, + int volume_number, int slice_number, + Ecat_field_name field, + int index, + int *ivalue, double *fvalue, char *svalue); +private int ecat_read_subhdr(Ecat_file *file, int volume, int slice); +private int ecat_lookup_field(Ecat_header_table_type *table, + Ecat_field_name field, + int *offset, int *length, Ecat_type *type, + char **description); +private void ecat_initialize_table(Ecat_header_table_type *table); +private int ecat_table_entry_compare(const void *v1, const void *v2); +private int ecat_table_offset_compare(const void *v1, const void *v2); +private int ecat_read_directory(Ecat_file *file); +private long get_dirblock(Ecat_file *file, long *dirblock, int offset); +private int ecat_get_subhdr_offset(Ecat_file *file, int volume, int slice, + long *offset); +private void *malloc_check(size_t size); +private void *realloc_check(void *ptr, size_t size); + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_open +@INPUT : filename - name of file to open +@OUTPUT : (none) +@RETURNS : Pointer to ECAT file descriptor or NULL if an error occurs. +@DESCRIPTION: Routine to open an ECAT file (for reading only), given + its pathname. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public Ecat_file *ecat_open(char *filename) +{ + Ecat_file *file; + + /* Allocate space for an ecat file structure */ + file = (void *) MALLOC(sizeof(*file)); + file->main_header = (void *) MALLOC(MAIN_HEADER_SIZE); + file->subheader = (void *) MALLOC(SUBHEADER_SIZE); + file->cur_subhdr_offset = -1; + file->subhdr_offsets = NULL; + + /* Open the file */ + if ((file->file_pointer=fopen(filename, "rb")) == NULL) { + ecat_close(file); + return NULL; + } + + /* Read in the main header */ + if (fread(file->main_header, sizeof(char), (size_t) MAIN_HEADER_SIZE, + file->file_pointer) != MAIN_HEADER_SIZE) { + ecat_close(file); + return NULL; + } + + /* Figure out which type of file we are using */ + if (strncmp((char *)file->main_header, MAGIC_STRING, + strlen(MAGIC_STRING)) == 0) { + file->header_description = ECAT_VER_7; + file->vax_byte_order = FALSE; + } + else { + file->header_description = ECAT_VER_PRE7; + file->vax_byte_order = TRUE; + } + + /* Get the number of frames, slices, bed positions and gates */ + if (ecat_get_main_value(file, ECAT_Num_Planes, 0, + &file->num_planes, NULL, NULL) || + ecat_get_main_value(file, ECAT_Num_Frames, 0, + &file->num_frames, NULL, NULL) || + ecat_get_main_value(file, ECAT_Num_Bed_Pos, 0, + &file->num_bed_positions, NULL, NULL) || + ecat_get_main_value(file, ECAT_Num_Gates, 0, + &file->num_gates, NULL, NULL)) { + ecat_close(file); + return NULL; + } + file->num_volumes = file->num_frames; + if (file->num_volumes < file->num_bed_positions) + file->num_volumes = file->num_bed_positions; + if (file->num_volumes < file->num_gates) + file->num_volumes = file->num_gates; + + /* Read the directory structure */ + if (ecat_read_directory(file)) { + ecat_close(file); + return NULL; + } + + /* Return the file pointer */ + return file; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_close +@INPUT : file - ecat file pointer +@OUTPUT : (none) +@RETURNS : (nothing) +@DESCRIPTION: Routine to close an ECAT file and free the associated structures +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public void ecat_close(Ecat_file *file) +{ + if (file==NULL) return; + if (file->file_pointer != NULL) { + (void) fclose(file->file_pointer); + } + if (file->subhdr_offsets != NULL) { + FREE(file->subhdr_offsets); + } + if (file->main_header != NULL) + FREE(file->main_header); + if (file->subheader != NULL) + FREE(file->subheader); + FREE(file); + + return; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_num_{planes,frames,bed_positions,gates} +@INPUT : file - ecat file pointer +@OUTPUT : (none) +@RETURNS : Number of * for an ECAT file. +@DESCRIPTION: Routine to get the number of planes (slices), frames, + bed_positions or gates in an ECAT file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_num_planes(Ecat_file *file) +{ + return file->num_planes; +} + +public int ecat_get_num_frames(Ecat_file *file) +{ + return file->num_frames; +} + +public int ecat_get_num_bed_positions(Ecat_file *file) +{ + return file->num_bed_positions; +} + +public int ecat_get_num_gates(Ecat_file *file) +{ + return file->num_gates; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_list_main +@INPUT : file - ecat file pointer + index - index into field list +@OUTPUT : (none) +@RETURNS : Next Ecat_field_name value, ECAT_No_Field if index is too large. +@DESCRIPTION: Routine to list the fields in an ECAT main header. + Should be called repeatedly with increasing values of index + (starting from 0) until field ECAT_No_Field is returned. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public Ecat_field_name ecat_list_main(Ecat_file *file, int index) +{ + return ecat_list_fields(file, ECAT_MAIN_HEADER, index); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_list_subhdr +@INPUT : file - ecat file pointer + index - index into field list +@OUTPUT : (none) +@RETURNS : Next Ecat_field_name value, ECAT_No_Field if index is too large. +@DESCRIPTION: Routine to list the fields in an ECAT subheader. + Should be called repeatedly with increasing values of index + (starting from 0) until field ECAT_No_Field is returned. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public Ecat_field_name ecat_list_subhdr(Ecat_file *file, int index) +{ + return ecat_list_fields(file, ECAT_SUBHEADER, index); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_list_fields +@INPUT : file - ecat file pointer + which_header - header in which to look + index - index into field list +@OUTPUT : (none) +@RETURNS : Next Ecat_field_name value, ECAT_No_Field if index is too large. +@DESCRIPTION: Routine to list the fields in an ECAT subheader. + Should be called repeatedly with increasing values of index + (starting from 0) until field ECAT_No_Field is returned. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private Ecat_field_name ecat_list_fields(Ecat_file *file, + Ecat_which_header which_header, + int index) +{ + Ecat_header_table_type *table; + + /* Get the appropriate header */ + switch (which_header) { + case ECAT_MAIN_HEADER: + table = file->header_description->main_header; + break; + case ECAT_SUBHEADER: + table = file->header_description->subheader; + break; + default: + return ECAT_No_Field; + } + + /* Initialize the table, if needed */ + if (!table->initialized) { + ecat_initialize_table(table); + } + + /* Check the index */ + if ((index < 0) || (index >= table->num_entries)) { + return ECAT_No_Field; + } + + /* Return the fields in file order */ + return table->file_order[index]; + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_main_field_length +@INPUT : file - ecat file pointer + field - ecat field name +@OUTPUT : (none) +@RETURNS : Length of field. Returns -1 if the field is not found. +@DESCRIPTION: Routine to get the length of a field from the ECAT main + header. Length is the number of elements, not the number of + bytes. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_main_field_length(Ecat_file *file, + Ecat_field_name field) +{ + return ecat_get_field_length(file, ECAT_MAIN_HEADER, field); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_subhdr_field_length +@INPUT : file - ecat file pointer + field - ecat field name +@OUTPUT : (none) +@RETURNS : Length of field. Returns -1 if the field is not found. +@DESCRIPTION: Routine to get the length of a field from an ECAT subheader. + Length is the number of elements, not the number of bytes. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_subhdr_field_length(Ecat_file *file, + Ecat_field_name field) +{ + return ecat_get_field_length(file, ECAT_SUBHEADER, field); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_field_length +@INPUT : file - ecat file pointer + which_header - header in which to look + field - ecat field name +@OUTPUT : (none) +@RETURNS : Length of field. Returns -1 if the field is not found. +@DESCRIPTION: Routine to get the length of a field from an ECAT + header. Length is the number of elements, not the number of + bytes. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_get_field_length(Ecat_file *file, + Ecat_which_header which_header, + Ecat_field_name field) +{ + Ecat_header_table_type *table; + Ecat_type type; + int length; + + /* Get the appropriate header */ + switch (which_header) { + case ECAT_MAIN_HEADER: + table = file->header_description->main_header; + break; + case ECAT_SUBHEADER: + table = file->header_description->subheader; + break; + default: + return -1; + } + + /* Look for the field description */ + if (ecat_lookup_field(table, field, NULL, &length, &type, NULL)) { + return -1; + } + + /* Return the length. If we have a string, then return 1. */ + if (type == ecat_char) + return 1; + else + return length; + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_main_field_description +@INPUT : file - ecat file pointer + field - ecat field name +@OUTPUT : (none) +@RETURNS : Description of field. Returns NULL if the field is not found. +@DESCRIPTION: Routine to get the description of a field from the ECAT main + header. The string returned should not be modified. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public char *ecat_get_main_field_description(Ecat_file *file, + Ecat_field_name field) +{ + return ecat_get_field_description(file, ECAT_MAIN_HEADER, field); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_subhdr_field_description +@INPUT : file - ecat file pointer + field - ecat field name +@OUTPUT : (none) +@RETURNS : Description of field. Returns NULL if the field is not found. +@DESCRIPTION: Routine to get the description of a field from an ECAT subheader. + The string returned should not be modified. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public char *ecat_get_subhdr_field_description(Ecat_file *file, + Ecat_field_name field) +{ + return ecat_get_field_description(file, ECAT_SUBHEADER, field); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_field_description +@INPUT : file - ecat file pointer + which_header - header in which to look + field - ecat field name +@OUTPUT : (none) +@RETURNS : Description of field. Returns NULL if the field is not found. +@DESCRIPTION: Routine to get the descrition of a field from an ECAT header. + The string returned should not be modified. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private char *ecat_get_field_description(Ecat_file *file, + Ecat_which_header which_header, + Ecat_field_name field) +{ + Ecat_header_table_type *table; + char *description; + + /* Get the appropriate header */ + switch (which_header) { + case ECAT_MAIN_HEADER: + table = file->header_description->main_header; + break; + case ECAT_SUBHEADER: + table = file->header_description->subheader; + break; + default: + return NULL; + } + + /* Look for the field description */ + if (ecat_lookup_field(table, field, NULL, NULL, NULL, &description)) { + return NULL; + } + + /* Return the description */ + return description; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_main_value +@INPUT : file - ecat file pointer + field - ecat field name + index - index for multi-valued fields (counting from zero) +@OUTPUT : ivalue - integer value + fvalue - floating-point value + svalue - string value +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to get a field value from the ECAT main header. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_main_value(Ecat_file *file, + Ecat_field_name field, int index, + int *ivalue, double *fvalue, char *svalue) +{ + return ecat_get_value(file, ECAT_MAIN_HEADER, 0, 0, field, index, + ivalue, fvalue, svalue); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_subhdr_value +@INPUT : file - ecat file pointer + volume - frame or bed position for subheader (from 0) + slice - slice number (counting from 0) + field - ecat field name + index - index for multi-valued fields (counting from zero) +@OUTPUT : ivalue - integer value + fvalue - floating-point value + svalue - string value +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to get a field value from the ECAT main header. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_subhdr_value(Ecat_file *file, int volume, int slice, + Ecat_field_name field, int index, + int *ivalue, double *fvalue, char *svalue) +{ + return ecat_get_value(file, ECAT_SUBHEADER, volume, slice, field, index, + ivalue, fvalue, svalue); +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_value +@INPUT : file - ECAT file pointer + which_header - header in which to look + volume_number - frame or bed position for subheader + (counting from 0) + slice_number - number of slice for subheader (counting from 0) + field - field to look for + index - index for multi-valued fields (counting from zero) +@OUTPUT : ivalue - integer value + fvalue - floating-point value + svalue - string value +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to look up a field value in a header table. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_get_value(Ecat_file *file, + Ecat_which_header which_header, + int volume_number, int slice_number, + Ecat_field_name field, + int index, + int *ivalue, double *fvalue, char *svalue) +{ + /* Sizes of ECAT types */ + static int ecat_type_sizes[] = {1, 2, 4, 4, 1}; + + /* Variables */ + unsigned char *header; + Ecat_header_table_type *table; + int offset, length; + Ecat_type type; + unsigned char byte_value; + short short_value; + long long_value; + float float_value; + char string[ECAT_MAX_STRING_LENGTH]; + + /* Get the appropriate header */ + switch (which_header) { + case ECAT_MAIN_HEADER: + header = file->main_header; + table = file->header_description->main_header; + break; + case ECAT_SUBHEADER: + if (ecat_read_subhdr(file, volume_number, slice_number)) { + return TRUE; + } + header = file->subheader; + table = file->header_description->subheader; + break; + default: + return TRUE; + } + + /* Look for the field description */ + if (ecat_lookup_field(table, field, &offset, &length, &type, NULL)) { + return TRUE; + } + + /* Check the index */ + if ((index < 0) || (index >= length) || + ((type == ecat_char) && (index > 0))) { + return TRUE; + } + offset += index * ecat_type_sizes[type]; + + /* Get the value and convert it */ + switch (type) { + case ecat_byte: + byte_value = header[offset]; + if (ivalue != NULL) *ivalue = byte_value; + if (fvalue != NULL) *fvalue = byte_value; + if (svalue != NULL) (void) sprintf(svalue, "%d", (int) byte_value); + break; + case ecat_short: + if (file->vax_byte_order) { + get_vax_short(1, &header[offset], &short_value); + } + else { + (void) memcpy(&short_value, &header[offset], sizeof(short)); + } + if (ivalue != NULL) *ivalue = short_value; + if (fvalue != NULL) *fvalue = short_value; + if (svalue != NULL) (void) sprintf(svalue, "%d", (int) short_value); + break; + case ecat_long: + if (file->vax_byte_order) { + get_vax_long(1, &header[offset], &long_value); + } + else { + (void) memcpy(&long_value, &header[offset], sizeof(long)); + } + if (ivalue != NULL) *ivalue = long_value; + if (fvalue != NULL) *fvalue = long_value; + if (svalue != NULL) (void) sprintf(svalue, "%d", (int) long_value); + break; + case ecat_float: + offset += index * sizeof(float); + if (file->vax_byte_order) { + get_vax_float(1, &header[offset], &float_value); + } + else { + (void) memcpy(&float_value, &header[offset], sizeof(float)); + } + if (ivalue != NULL) *ivalue = float_value; + if (fvalue != NULL) *fvalue = float_value; + if (svalue != NULL) (void) sprintf(svalue, "%.7g", + (double) float_value); + break; + case ecat_char: + (void) memcpy(string, &header[offset], length); + string[length] = '\0'; + if (ivalue != NULL) *ivalue = atoi(svalue); + if (fvalue != NULL) *fvalue = atof(svalue); + if (svalue != NULL) (void) strcpy(svalue, string); + break; + default: + return TRUE; + } + + return FALSE; + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_read_subhdr +@INPUT : file - ecat file pointer + volume - frame or bed position for subheader (from 0) + slice - slice number (counting from 0) +@OUTPUT : file - file structure is modified to include new subheader +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to read in the appropriate subheader, if needed. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_read_subhdr(Ecat_file *file, int volume, int slice) +{ + long offset; + + /* Get subheader offset */ + if (ecat_get_subhdr_offset(file, volume, slice, &offset)) { + return TRUE; + } + + /* Do we need to read in the header? */ + if (offset == file->cur_subhdr_offset) { + return FALSE; + } + + /* If so, do it */ + if (fseek(file->file_pointer, offset, SEEK_SET) || + (fread(file->subheader, sizeof(char), (size_t) SUBHEADER_SIZE, + file->file_pointer) != SUBHEADER_SIZE)) { + return TRUE; + } + file->cur_subhdr_offset = offset; + + return FALSE; + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_lookup_field +@INPUT : table - table containing header fields + field - field to look for +@OUTPUT : offset - offset into header + length - length of field (number of values, not bytes) + type - type of values + description - string describing field (ptr is passed back + by reference) +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to look up a field value in a header table. All + output arguments can be NULL to indicate that the + corresponding value should not be returned. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_lookup_field(Ecat_header_table_type *table, + Ecat_field_name field, + int *offset, int *length, Ecat_type *type, + char **description) +{ + Ecat_field_description_type field_description, *result; + + /* Initialize the table, if needed */ + if (!table->initialized) { + ecat_initialize_table(table); + } + + /* Set up a dummy entry */ + field_description.name = field; + + /* Search for the field */ + result = bsearch(&field_description, table->fields, table->num_entries, + sizeof(table->fields[0]), ecat_table_entry_compare); + + if (result == NULL) { + return TRUE; + } + + /* Store the results */ + if (offset != NULL) *offset = result->offset; + if (length != NULL) *length = result->length; + if (type != NULL) *type = result->type; + if (description != NULL) *description = result->description; + + return FALSE; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_initialize_table +@INPUT : table - table containing header fields +@OUTPUT : table - modified appropriately +@RETURNS : (nothing) +@DESCRIPTION: Routine to initialize the header table. +@METHOD : Initialized should be set to TRUE, num_entries should already + be set and fields should point to a list. The fields list will + be sorted and the file_order field will be set. +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private void ecat_initialize_table(Ecat_header_table_type *table) +{ + Ecat_field_description_type *file_order_list; + int iel; + + /* Check that table is not already initialized */ + if (table->initialized) return; + + /* Sort the field list by field name for binary searching */ + qsort(table->fields, table->num_entries, sizeof(table->fields[0]), + ecat_table_entry_compare); + + /* Get a list of field names in file order by copying the field list + and sorting it by offset */ + file_order_list = MALLOC(table->num_entries * sizeof(*file_order_list)); + for (iel=0; iel < table->num_entries; iel++) { + file_order_list[iel] = table->fields[iel]; + } + qsort(file_order_list, table->num_entries, sizeof(file_order_list[0]), + ecat_table_offset_compare); + table->file_order = MALLOC(table->num_entries * + sizeof(table->file_order[0])); + for (iel=0; iel < table->num_entries; iel++) { + table->file_order[iel] = file_order_list[iel].name; + } + FREE(file_order_list); + + /* Mark the table as initialized */ + table->initialized = TRUE; + + return; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_table_entry_compare +@INPUT : v1 - first value +@OUTPUT : v2 - second value +@RETURNS : (nothing) +@DESCRIPTION: Routine to compare field list elements for sorting by name. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_table_entry_compare(const void *v1, const void *v2) +{ + Ecat_field_description_type *first, *second; + + /* Get pointers */ + first = (Ecat_field_description_type *) v1; + second = (Ecat_field_description_type *) v2; + + /* Compare field names */ + return ((int) first->name - (int) second->name); + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_table_offset_compare +@INPUT : v1 - first value +@OUTPUT : v2 - second value +@RETURNS : (nothing) +@DESCRIPTION: Routine to compare field list elements for sorting by offset. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_table_offset_compare(const void *v1, const void *v2) +{ + Ecat_field_description_type *first, *second; + + /* Get pointers */ + first = (Ecat_field_description_type *) v1; + second = (Ecat_field_description_type *) v2; + + /* Compare field names */ + return ((int) first->offset - (int) second->offset); + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_read_directory +@INPUT : file - ecat file pointer +@OUTPUT : (nothing) +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to read in the ECAT subheader directory +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_read_directory(Ecat_file *file) +{ + int ientry, start_entry, num_used, num_alloc; + long nextblock, dirblock[DIRBLOCK_SIZE / sizeof(long)]; + + /* Allocate space for the subheader offset array */ + num_alloc = file->num_volumes; + if (file->header_description == ECAT_VER_PRE7) + num_alloc *= file->num_planes; + file->num_subhdrs = 0; + file->subhdr_offsets = + MALLOC(num_alloc * sizeof(file->subhdr_offsets[0])); + + /* Reading directory blocks until done */ + nextblock = FIRST_DIRBLOCK; + do { + + /* Read in the block */ + if (fseek(file->file_pointer, (nextblock - 1) * BLOCK_SIZE, SEEK_SET) || + (fread(dirblock, sizeof(char), sizeof(dirblock), + file->file_pointer) != sizeof(dirblock))) { + return TRUE; + } + + /* Get a pointer to the next block and the number of entries used */ + nextblock = get_dirblock(file, dirblock, DRBLK_NEXT); + num_used = get_dirblock(file, dirblock, DRBLK_NUMUSED); + + /* Increment the number of subheaders */ + start_entry = file->num_subhdrs; + file->num_subhdrs += num_used; + if (num_alloc < file->num_subhdrs) { + num_alloc = file->num_subhdrs; + REALLOC(file->subhdr_offsets, + num_alloc * sizeof(file->subhdr_offsets[0])); + } + + /* Save the offsets */ + for (ientry=0; ientry < num_used; ientry++) { + file->subhdr_offsets[start_entry + ientry] = BLOCK_SIZE * + (get_dirblock(file, dirblock, + DRBLK_SKIP + ientry * DRBLK_WIDTH + DRBLK_PTR) - 1); + } + + } while (nextblock > FIRST_DIRBLOCK); + + return FALSE; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : get_dirblock +@INPUT : file - ecat file pointer + dirblock - directory block + offset - offset (in longwords) into the block +@OUTPUT : (nothing) +@RETURNS : directory block value +@DESCRIPTION: Routine to get a value from an ECAT directory block. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private long get_dirblock(Ecat_file *file, long *dirblock, int offset) +{ + long value; + + if (file->header_description == ECAT_VER_PRE7) { + get_vax_long(1, &dirblock[offset], &value); + } + else if (file->header_description == ECAT_VER_7) { + value = dirblock[offset]; + } + else { + return 0; + } + + return value; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_image +@INPUT : file - ecat file pointer + volume - frame or bed position (from 0) + slice - slice number (counting from 0) +@OUTPUT : image +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to get an image from an ECAT file +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +public int ecat_get_image(Ecat_file *file, int volume, int slice, + short *image) +{ + long file_offset; + int xsize, ysize, zsize, data_type, bytes_per_pixel; + long image_npix, image_size, array_offset, ipix; + unsigned char *bimage; + + /* Get the image size and type */ + if (ecat_get_subhdr_value(file, volume, slice, ECAT_X_Dimension, 0, + &xsize, NULL, NULL) || + ecat_get_subhdr_value(file, volume, slice, ECAT_Y_Dimension, 0, + &ysize, NULL, NULL) || + ecat_get_subhdr_value(file, volume, slice, ECAT_Data_Type, 0, + &data_type, NULL, NULL) || + (xsize <= 0) || (ysize <= 0)) { + return TRUE; + } + + /* Figure out the image size */ + switch (data_type) { + case 1: + bytes_per_pixel = 1; break; + case 2: + case 6: + bytes_per_pixel = 2; break; + default: + return TRUE; + } + image_npix = xsize * ysize; + image_size = image_npix * bytes_per_pixel; + + /* Look for a z size */ + zsize = 0; + (void) ecat_get_subhdr_value(file, volume, slice, ECAT_Z_Dimension, 0, + &zsize, NULL, NULL); + + /* Check the that the slice is in range */ + if ((slice < 0) || ((zsize > 0) && (slice > zsize))) { + return TRUE; + } + + /* Find the appropriate subheader */ + if (ecat_get_subhdr_offset(file, volume, slice, &file_offset)) { + return TRUE; + } + + /* Adjust the offset appropriately */ + file_offset += BLOCK_SIZE; + if (zsize > 0) { + file_offset += image_size * slice; + } + + /* Calculate image size and offsets */ + array_offset = image_npix * (sizeof(short) - bytes_per_pixel); + bimage = (unsigned char *) image; + + /* Read in the image */ + if (fseek(file->file_pointer, file_offset, SEEK_SET) || + (fread(&bimage[array_offset], (size_t) bytes_per_pixel, + (size_t) image_npix, file->file_pointer) != image_npix)) { + return TRUE; + } + + /* Transform the image to the right type */ + switch (bytes_per_pixel) { + case 1: + for (ipix=0; ipix<image_npix; ipix++) { + image[ipix] = bimage[array_offset+ipix]; + } + break; + case 2: + if (file->header_description == ECAT_VER_PRE7) { + get_vax_short(image_npix, image, image); + } + break; + default: + return TRUE; + } + + return FALSE; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_get_subhdr_offset +@INPUT : file - ecat file pointer + volume - frame or bed position for subheader (from 0) + slice - slice number (counting from 0) +@OUTPUT : offset - offset in bytes to appropriate subheader +@RETURNS : FALSE if successful, TRUE otherwise +@DESCRIPTION: Routine to calculate the offset to the appropriate subheader. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private int ecat_get_subhdr_offset(Ecat_file *file, int volume, int slice, + long *offset) +{ + int subheader_number; + + /* Check to see if volume and slice are within bounds */ + if ((volume < 0) || (volume >= file->num_volumes) || + (slice < 0) || (slice >= file->num_planes)) { + return TRUE; + } + + /* Calculate subheader number */ + if (file->header_description == ECAT_VER_7) { + subheader_number = volume; + } + else if (file->header_description == ECAT_VER_PRE7) { + subheader_number = volume * file->num_planes + slice; + } + else { + return TRUE; + } + + /* Check the subheader number */ + if ((subheader_number < 0) || (subheader_number >= file->num_subhdrs)) { + return TRUE; + } + + /* Get the offset for this subheader */ + *offset = file->subhdr_offsets[subheader_number]; + + return FALSE; + +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : malloc_check +@INPUT : size - number of bytes to allocate +@OUTPUT : (nothing) +@RETURNS : pointer to memory +@DESCRIPTION: Routine to allocate memory. It will never return NULL - the + program will exit first. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private void *malloc_check(size_t size) +{ + void *ptr; + + ptr = malloc(size); + if (ptr == NULL) { + (void) fprintf(stderr, "Out of memory.\n"); + exit(EXIT_FAILURE); + } + return ptr; +} + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : realloc_check +@INPUT : ptr - old pointer + size - number of bytes to allocate +@OUTPUT : (nothing) +@RETURNS : pointer to memory +@DESCRIPTION: Routine to re-allocate memory. It will never return NULL - the + program will exit first. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +private void *realloc_check(void *ptr, size_t size) +{ + ptr = realloc(ptr, size); + if (ptr == NULL) { + (void) fprintf(stderr, "Out of memory.\n"); + exit(EXIT_FAILURE); + } + return ptr; +} +
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/ecat_file.h @@ -0,0 +1,172 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_file.h +@DESCRIPTION: Header file for routines that read ECAT image files +@GLOBALS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : $Log: ecat_file.h,v $ +@MODIFIED : Revision 1.1 1996-01-18 14:52:14 neelin +@MODIFIED : Initial revision +@MODIFIED : +@COPYRIGHT : + Copyright 1996 Peter Neelin, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ + +#ifndef public +#define public +#endif +#ifndef private +#define private static +#endif + +#define ECAT_MAX_STRING_LENGTH 64 + +typedef enum { + ECAT_No_Field, + ECAT_Magic_Number, + ECAT_Original_Filename, + ECAT_Sw_Version, + ECAT_System_Type, + ECAT_File_Type, + ECAT_Serial_Number, + ECAT_Scan_Start_Time, + ECAT_Isotope_Name, + ECAT_Isotope_Halflife, + ECAT_Radiopharmaceutical, + ECAT_Gantry_Tilt, + ECAT_Gantry_Rotation, + ECAT_Bed_Elevation, + ECAT_Intrinsic_Tilt, + ECAT_Wobble_Speed, + ECAT_Transm_Source_Type, + ECAT_Distance_Scanned, + ECAT_Transaxial_Fov, + ECAT_Angular_Compression, + ECAT_Coin_Samp_Mode, + ECAT_Axial_Samp_Mode, + ECAT_Calibration_Factor, + ECAT_Calibration_Units, + ECAT_Calibration_Units_Label, + ECAT_Compression_Code, + ECAT_Study_Type, + ECAT_Patient_Id, + ECAT_Patient_Name, + ECAT_Patient_Sex, + ECAT_Patient_Dexterity, + ECAT_Patient_Age, + ECAT_Patient_Height, + ECAT_Patient_Weight, + ECAT_Patient_Birth_Date, + ECAT_Physician_Name, + ECAT_Operator_Name, + ECAT_Study_Description, + ECAT_Acquision_Type, + ECAT_Patient_Orientation, + ECAT_Facility_Name, + ECAT_Num_Planes, + ECAT_Num_Frames, + ECAT_Num_Gates, + ECAT_Num_Bed_Pos, + ECAT_Init_Bed_Position, + ECAT_Bed_Position, + ECAT_Plane_Separation, + ECAT_Lwr_Sctr_Thres, + ECAT_Lwr_True_Thres, + ECAT_Upr_True_Thres, + ECAT_User_Process_Code, + ECAT_Acquisition_Mode, + ECAT_Bin_Size, + ECAT_Branching_Fraction, + ECAT_Dose_Start_Time, + ECAT_Dosage, + ECAT_Well_Counter_Corr_Factor, + ECAT_Data_Units, + ECAT_Septa_State, + ECAT_Data_Type, + ECAT_Num_Dimensions, + ECAT_X_Dimension, + ECAT_Y_Dimension, + ECAT_Z_Dimension, + ECAT_X_Offset, + ECAT_Y_Offset, + ECAT_Z_Offset, + ECAT_Recon_Zoom, + ECAT_Scale_Factor, + ECAT_Image_Min, + ECAT_Image_Max, + ECAT_X_Pixel_Size, + ECAT_Y_Pixel_Size, + ECAT_Z_Pixel_Size, + ECAT_Frame_Duration, + ECAT_Frame_Start_Time, + ECAT_Filter_Code, + ECAT_X_Resolution, + ECAT_Y_Resolution, + ECAT_Z_Resolution, + ECAT_X_Rotation_Angle, + ECAT_Y_Rotation_Angle, + ECAT_Z_Rotation_Angle, + ECAT_Decay_Corr_Fctr, + ECAT_Corrections_Applied, + ECAT_Gate_Duration, + ECAT_R_Wave_Offset, + ECAT_Num_Accepted_Beats, + ECAT_Filter_Cutoff_Frequence, + ECAT_Filter_Dc_Component, + ECAT_Filter_Ramp_Slope, + ECAT_Filter_Order, + ECAT_Filter_Scatter_Fraction, + ECAT_Filetr_Scatter_Slope, + ECAT_Annotation, + ECAT_Da_X_Rotation_Angle, + ECAT_Da_Y_Rotation_Angle, + ECAT_Da_Z_Rotation_Angle, + ECAT_Da_X_Translation, + ECAT_Da_Y_Translation, + ECAT_Da_Z_Translation, + ECAT_Da_X_Scale_Factor, + ECAT_Da_Y_Scale_Factor, + ECAT_Da_Z_Scale_Factor, + ECAT_Scan_Start_Day, + ECAT_Scan_Start_Month, + ECAT_Scan_Start_Year, + ECAT_Scan_Start_Hour, + ECAT_Scan_Start_Minute, + ECAT_Scan_Start_Second, + ECAT_Rot_Source_Speed +} Ecat_field_name; + +typedef struct Ecat_file Ecat_file; + +/* Routine declarations */ +public Ecat_file *ecat_open(char *filename); +public void ecat_close(Ecat_file *file); +public int ecat_get_num_planes(Ecat_file *file); +public int ecat_get_num_frames(Ecat_file *file); +public int ecat_get_num_bed_positions(Ecat_file *file); +public int ecat_get_num_gates(Ecat_file *file); +public Ecat_field_name ecat_list_main(Ecat_file *file, int index); +public Ecat_field_name ecat_list_subhdr(Ecat_file *file, int index); +public int ecat_get_main_field_length(Ecat_file *file, + Ecat_field_name field); +public int ecat_get_subhdr_field_length(Ecat_file *file, + Ecat_field_name field); +public char *ecat_get_main_field_description(Ecat_file *file, + Ecat_field_name field); +public char *ecat_get_subhdr_field_description(Ecat_file *file, + Ecat_field_name field); +public int ecat_get_main_value(Ecat_file *file, + Ecat_field_name field, int index, + int *ivalue, double *fvalue, char *svalue); +public int ecat_get_subhdr_value(Ecat_file *file, int volume, int slice, + Ecat_field_name field, int index, + int *ivalue, double *fvalue, char *svalue); +public int ecat_get_image(Ecat_file *file, int volume, int slice, + short *image);
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/ecat_header_definition.h @@ -0,0 +1,244 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecat_header_definition.h +@DESCRIPTION: Header file containing header description for ECAT files +@GLOBALS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : $Log: ecat_header_definition.h,v $ +@MODIFIED : Revision 1.1 1996-01-18 14:52:14 neelin +@MODIFIED : Initial revision +@MODIFIED : +@COPYRIGHT : + Copyright 1996 Peter Neelin, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ + +/* Macro to set up header definition tables */ +#define INITIAL_HEADER_TABLE(field_list) \ +{FALSE, sizeof(field_list)/sizeof(field_list[0]), field_list, NULL} + +/* Main header and sub header descriptions for ECAT version 7 files */ +static Ecat_field_description_type version_7_main_list[] = +{ + {ECAT_Magic_Number, 0, 14, ecat_char, "Magic_Number"}, + {ECAT_Original_Filename, 14, 32, ecat_char, "Original_Filename"}, + {ECAT_Sw_Version, 46, 1, ecat_short, "Sw_Version"}, + {ECAT_System_Type, 48, 1, ecat_short, "System_Type"}, + {ECAT_File_Type, 50, 1, ecat_short, "File_Type"}, + {ECAT_Serial_Number, 52, 10, ecat_char, "Serial_Number"}, + {ECAT_Scan_Start_Time, 62, 1, ecat_long, "Scan_Start_Time"}, + {ECAT_Isotope_Name, 66, 8, ecat_char, "Isotope_Name"}, + {ECAT_Isotope_Halflife, 74, 1, ecat_float, "Isotope_Halflife"}, + {ECAT_Radiopharmaceutical, 78, 32, ecat_char, "Radiopharmaceutical"}, + {ECAT_Gantry_Tilt, 110, 1, ecat_float, "Gantry_Tilt"}, + {ECAT_Gantry_Rotation, 114, 1, ecat_float, "Gantry_Rotation"}, + {ECAT_Bed_Elevation, 118, 1, ecat_float, "Bed_Elevation"}, + {ECAT_Intrinsic_Tilt, 122, 1, ecat_float, "Intrinsic_Tilt"}, + {ECAT_Wobble_Speed, 126, 1, ecat_short, "Wobble_Speed"}, + {ECAT_Transm_Source_Type, 128, 1, ecat_short, "Transm_Source_Type"}, + {ECAT_Distance_Scanned, 130, 1, ecat_float, "Distance_Scanned"}, + {ECAT_Transaxial_Fov, 134, 1, ecat_float, "Transaxial_Fov"}, + {ECAT_Angular_Compression, 138, 1, ecat_short, "Angular_Compression"}, + {ECAT_Coin_Samp_Mode, 140, 1, ecat_short, "Coin_Samp_Mode"}, + {ECAT_Axial_Samp_Mode, 142, 1, ecat_short, "Axial_Samp_Mode"}, + {ECAT_Calibration_Factor, 144, 1, ecat_float, "Calibration_Factor"}, + {ECAT_Calibration_Units, 148, 1, ecat_short, "Calibration_Units"}, + {ECAT_Calibration_Units_Label, 150, 1, ecat_short, "Calibration_Units_Label"}, + {ECAT_Compression_Code, 152, 1, ecat_short, "Compression_Code"}, + {ECAT_Study_Type, 154, 12, ecat_char, "Study_Type"}, + {ECAT_Patient_Id, 166, 16, ecat_char, "Patient_Id"}, + {ECAT_Patient_Name, 182, 32, ecat_char, "Patient_Name"}, + {ECAT_Patient_Sex, 214, 1, ecat_char, "Patient_Sex"}, + {ECAT_Patient_Dexterity, 215, 1, ecat_char, "Patient_Dexterity"}, + {ECAT_Patient_Age, 216, 1, ecat_float, "Patient_Age"}, + {ECAT_Patient_Height, 220, 1, ecat_float, "Patient_Height"}, + {ECAT_Patient_Weight, 224, 1, ecat_float, "Patient_Weight"}, + {ECAT_Patient_Birth_Date, 228, 1, ecat_long, "Patient_Birth_Date"}, + {ECAT_Physician_Name, 232, 32, ecat_char, "Physician_Name"}, + {ECAT_Operator_Name, 264, 32, ecat_char, "Operator_Name"}, + {ECAT_Study_Description, 296, 32, ecat_char, "Study_Description"}, + {ECAT_Acquision_Type, 328, 1, ecat_short, "Acquision_Type"}, + {ECAT_Patient_Orientation, 330, 1, ecat_short, "Patient_Orientation"}, + {ECAT_Facility_Name, 332, 20, ecat_char, "Facility_Name"}, + {ECAT_Num_Planes, 352, 1, ecat_short, "Num_Planes"}, + {ECAT_Num_Frames, 354, 1, ecat_short, "Num_Frames"}, + {ECAT_Num_Gates, 356, 1, ecat_short, "Num_Gates"}, + {ECAT_Num_Bed_Pos, 358, 1, ecat_short, "Num_Bed_Pos"}, + {ECAT_Init_Bed_Position, 360, 1, ecat_float, "Init_Bed_Position"}, + {ECAT_Bed_Position, 364, 15, ecat_float, "Bed_Position"}, + {ECAT_Plane_Separation, 424, 1, ecat_float, "Plane_Separation"}, + {ECAT_Lwr_Sctr_Thres, 428, 1, ecat_short, "Lwr_Sctr_Thres"}, + {ECAT_Lwr_True_Thres, 430, 1, ecat_short, "Lwr_True_Thres"}, + {ECAT_Upr_True_Thres, 432, 1, ecat_short, "Upr_True_Thres"}, + {ECAT_User_Process_Code, 434, 10, ecat_char, "User_Process_Code"}, + {ECAT_Acquisition_Mode, 444, 1, ecat_short, "Acquisition_Mode"}, + {ECAT_Bin_Size, 446, 1, ecat_float, "Bin_Size"}, + {ECAT_Branching_Fraction, 450, 1, ecat_float, "Branching_Fraction"}, + {ECAT_Dose_Start_Time, 454, 1, ecat_long, "Dose_Start_Time"}, + {ECAT_Dosage, 458, 1, ecat_float, "Dosage"}, + {ECAT_Well_Counter_Corr_Factor, 462, 1, ecat_float, "Well_Counter_Corr_Factor"} +}; + +static Ecat_field_description_type version_7_subhdr_list[] = +{ + {ECAT_Data_Type, 0, 1, ecat_short, "Data_Type"}, + {ECAT_Num_Dimensions, 2, 1, ecat_short, "Num_Dimensions"}, + {ECAT_X_Dimension, 4, 1, ecat_short, "X_Dimension"}, + {ECAT_Y_Dimension, 6, 1, ecat_short, "Y_Dimension"}, + {ECAT_Z_Dimension, 8, 1, ecat_short, "Z_Dimension"}, + {ECAT_X_Offset, 10, 1, ecat_float, "X_Offset"}, + {ECAT_Y_Offset, 14, 1, ecat_float, "Y_Offset"}, + {ECAT_Z_Offset, 18, 1, ecat_float, "Z_Offset"}, + {ECAT_Recon_Zoom, 22, 1, ecat_float, "Recon_Zoom"}, + {ECAT_Scale_Factor, 26, 1, ecat_float, "Scale_Factor"}, + {ECAT_Image_Min, 30, 1, ecat_short, "Image_Min"}, + {ECAT_Image_Max, 32, 1, ecat_short, "Image_Max"}, + {ECAT_X_Pixel_Size, 34, 1, ecat_float, "X_Pixel_Size"}, + {ECAT_Y_Pixel_Size, 38, 1, ecat_float, "Y_Pixel_Size"}, + {ECAT_Z_Pixel_Size, 42, 1, ecat_float, "Z_Pixel_Size"}, + {ECAT_Frame_Duration, 46, 1, ecat_long, "Frame_Duration"}, + {ECAT_Frame_Start_Time, 50, 1, ecat_long, "Frame_Start_Time"}, + {ECAT_Filter_Code, 54, 1, ecat_short, "Filter_Code"}, + {ECAT_X_Resolution, 56, 1, ecat_float, "X_Resolution"}, + {ECAT_Y_Resolution, 60, 1, ecat_float, "Y_Resolution"}, + {ECAT_Z_Resolution, 64, 1, ecat_float, "Z_Resolution"}, + {ECAT_X_Rotation_Angle, 68, 1, ecat_float, "X_Rotation_Angle"}, + {ECAT_Y_Rotation_Angle, 72, 1, ecat_float, "Y_Rotation_Angle"}, + {ECAT_Z_Rotation_Angle, 76, 1, ecat_float, "Z_Rotation_Angle"}, + {ECAT_Decay_Corr_Fctr, 80, 1, ecat_float, "Decay_Corr_Fctr"}, + {ECAT_Corrections_Applied, 84, 1, ecat_long, "Corrections_Applied"}, + {ECAT_Gate_Duration, 88, 1, ecat_long, "Gate_Duration"}, + {ECAT_R_Wave_Offset, 92, 1, ecat_long, "R_Wave_Offset"}, + {ECAT_Num_Accepted_Beats, 96, 1, ecat_long, "Num_Accepted_Beats"}, + {ECAT_Filter_Cutoff_Frequence, 100, 1, ecat_float, "Filter_Cutoff_Frequence"}, + {ECAT_Filter_Dc_Component, 104, 1, ecat_float, "Filter_Dc_Component"}, + {ECAT_Filter_Ramp_Slope, 108, 1, ecat_float, "Filter_Ramp_Slope"}, + {ECAT_Filter_Order, 112, 1, ecat_short, "Filter_Order"}, + {ECAT_Filter_Scatter_Fraction, 114, 1, ecat_float, "Filter_Scatter_Fraction"}, + {ECAT_Filetr_Scatter_Slope, 118, 1, ecat_float, "Filetr_Scatter_Slope"}, + {ECAT_Annotation, 122, 40, ecat_char, "Annotation"}, + {ECAT_Da_X_Rotation_Angle, 162, 1, ecat_float, "Da_X_Rotation_Angle"}, + {ECAT_Da_Y_Rotation_Angle, 166, 1, ecat_float, "Da_Y_Rotation_Angle"}, + {ECAT_Da_Z_Rotation_Angle, 170, 1, ecat_float, "Da_Z_Rotation_Angle"}, + {ECAT_Da_X_Translation, 174, 1, ecat_float, "Da_X_Translation"}, + {ECAT_Da_Y_Translation, 178, 1, ecat_float, "Da_Y_Translation"}, + {ECAT_Da_Z_Translation, 182, 1, ecat_float, "Da_Z_Translation"}, + {ECAT_Da_X_Scale_Factor, 186, 1, ecat_float, "Da_X_Scale_Factor"}, + {ECAT_Da_Y_Scale_Factor, 190, 1, ecat_float, "Da_Y_Scale_Factor"}, + {ECAT_Da_Z_Scale_Factor, 194, 1, ecat_float, "Da_Z_Scale_Factor"} +}; + +static Ecat_header_table_type +version_7_main_table = INITIAL_HEADER_TABLE(version_7_main_list); + +static Ecat_header_table_type +version_7_subhdr_table = INITIAL_HEADER_TABLE(version_7_subhdr_list); + +static Ecat_header_description_type version_7_description = +{ + &version_7_main_table, + &version_7_subhdr_table +}; + +static Ecat_header_description_type *ECAT_VER_7 = &version_7_description; + +/* Main header and sub header descriptions for ECAT files before version 7 */ +static Ecat_field_description_type version_pre7_main_list[] = +{ + {ECAT_Original_Filename, 28, 20, ecat_char, "Original_Filename"}, + {ECAT_Sw_Version, 48, 1, ecat_short, "Sw_Version"}, + {ECAT_Data_Type, 50, 1, ecat_short, "Data_type"}, + {ECAT_System_Type, 52, 1, ecat_short, "System_Type"}, + {ECAT_File_Type, 54, 1, ecat_short, "File_Type"}, + {ECAT_Serial_Number, 56, 10, ecat_char, "Serial_Number"}, + {ECAT_Scan_Start_Day, 66, 1, ecat_short, "Scan_Start_Day"}, + {ECAT_Scan_Start_Month, 68, 1, ecat_short, "Scan_Start_Month"}, + {ECAT_Scan_Start_Year, 70, 1, ecat_short, "Scan_Start_Year"}, + {ECAT_Scan_Start_Hour, 72, 1, ecat_short, "Scan_Start_Hour"}, + {ECAT_Scan_Start_Minute, 74, 1, ecat_short, "Scan_Start_Minute"}, + {ECAT_Scan_Start_Second, 76, 1, ecat_short, "Scan_Start_Second"}, + {ECAT_Isotope_Name, 78, 8, ecat_char, "Isotope_Name"}, + {ECAT_Isotope_Halflife, 86, 1, ecat_float, "Isotope_Halflife"}, + {ECAT_Radiopharmaceutical, 90, 32, ecat_char, "Radiopharmaceutical"}, + {ECAT_Gantry_Tilt, 122, 1, ecat_float, "Gantry_Tilt"}, + {ECAT_Gantry_Rotation, 126, 1, ecat_float, "Gantry_Rotation"}, + {ECAT_Bed_Elevation, 130, 1, ecat_float, "Bed_Elevation"}, + {ECAT_Rot_Source_Speed, 134, 1, ecat_short, "Rot_Source_Speed"}, + {ECAT_Wobble_Speed, 136, 1, ecat_short, "Wobble_Speed"}, + {ECAT_Transm_Source_Type, 138, 1, ecat_short, "Transm_Source_Type"}, + {ECAT_Distance_Scanned, 140, 1, ecat_float, "Distance_Scanned"}, + {ECAT_Transaxial_Fov, 144, 1, ecat_float, "Transaxial_Fov"}, + {ECAT_Coin_Samp_Mode, 150, 1, ecat_short, "Coin_Samp_Mode"}, + {ECAT_Axial_Samp_Mode, 152, 1, ecat_short, "Axial_Samp_Mode"}, + {ECAT_Calibration_Factor, 154, 1, ecat_float, "Calibration_Factor"}, + {ECAT_Calibration_Units, 158, 1, ecat_short, "Calibration_Units"}, + {ECAT_Compression_Code, 160, 1, ecat_short, "Compression_Code"}, + {ECAT_Study_Type, 162, 12, ecat_char, "Study_Type"}, + {ECAT_Patient_Id, 174, 16, ecat_char, "Patient_Id"}, + {ECAT_Patient_Name, 190, 32, ecat_char, "Patient_Name"}, + {ECAT_Patient_Sex, 222, 1, ecat_char, "Patient_Sex"}, + {ECAT_Patient_Age, 223, 10, ecat_char, "Patient_Age"}, + {ECAT_Patient_Height, 233, 10, ecat_char, "Patient_Height"}, + {ECAT_Patient_Weight, 243, 10, ecat_char, "Patient_Weight"}, + {ECAT_Patient_Dexterity, 253, 1, ecat_char, "Patient_Dexterity"}, + {ECAT_Physician_Name, 254, 32, ecat_char, "Physician_Name"}, + {ECAT_Operator_Name, 286, 32, ecat_char, "Operator_Name"}, + {ECAT_Study_Description, 318, 32, ecat_char, "Study_Description"}, + {ECAT_Acquision_Type, 350, 1, ecat_short, "Acquision_Type"}, + {ECAT_Facility_Name, 356, 20, ecat_char, "Facility_Name"}, + {ECAT_Num_Planes, 376, 1, ecat_short, "Num_Planes"}, + {ECAT_Num_Frames, 378, 1, ecat_short, "Num_Frames"}, + {ECAT_Num_Gates, 380, 1, ecat_short, "Num_Gates"}, + {ECAT_Num_Bed_Pos, 382, 1, ecat_short, "Num_Bed_Pos"}, + {ECAT_Init_Bed_Position, 384, 1, ecat_float, "Init_Bed_Position"}, + {ECAT_Bed_Position, 388, 15, ecat_float, "Bed_Position"}, + {ECAT_Plane_Separation, 448, 1, ecat_float, "Plane_Separation"}, + {ECAT_Lwr_Sctr_Thres, 452, 1, ecat_short, "Lwr_Sctr_Thres"}, + {ECAT_Lwr_True_Thres, 454, 1, ecat_short, "Lwr_True_Thres"}, + {ECAT_Upr_True_Thres, 456, 1, ecat_short, "Upr_True_Thres"}, + {ECAT_User_Process_Code, 462, 10, ecat_char, "User_Process_Code"} +}; + +static Ecat_field_description_type version_pre7_subhdr_list[] = +{ + {ECAT_Data_Type, 126, 1, ecat_short, "Data_Type"}, + {ECAT_Num_Dimensions, 128, 1, ecat_short, "Num_Dimensions"}, + {ECAT_X_Dimension, 132, 1, ecat_short, "X_Dimension"}, + {ECAT_Y_Dimension, 134, 1, ecat_short, "Y_Dimension"}, + {ECAT_X_Offset, 160, 1, ecat_float, "X_Offset"}, + {ECAT_Y_Offset, 164, 1, ecat_float, "Y_Offset"}, + {ECAT_Recon_Zoom, 168, 1, ecat_float, "Recon_Zoom"}, + {ECAT_Scale_Factor, 172, 1, ecat_float, "Scale_Factor"}, + {ECAT_Image_Min, 176, 1, ecat_short, "Image_Min"}, + {ECAT_Image_Max, 178, 1, ecat_short, "Image_Max"}, + {ECAT_X_Pixel_Size, 184, 1, ecat_float, "X_Pixel_Size"}, + {ECAT_Y_Pixel_Size, 184, 1, ecat_float, "Y_Pixel_Size"}, /* Same location */ + {ECAT_Z_Pixel_Size, 188, 1, ecat_float, "Z_Pixel_Size"}, + {ECAT_Frame_Duration, 192, 1, ecat_long, "Frame_Duration"}, + {ECAT_Frame_Start_Time, 196, 1, ecat_long, "Frame_Start_Time"}, + {ECAT_Filter_Code, 236, 1, ecat_short, "Filter_Code"}, + {ECAT_Z_Rotation_Angle, 296, 1, ecat_float, "Z_Rotation_Angle"}, + {ECAT_Decay_Corr_Fctr, 304, 1, ecat_float, "Decay_Corr_Fctr"}, + {ECAT_Annotation, 420, 40, ecat_char, "Annotation"} +}; + +static Ecat_header_table_type +version_pre7_main_table = INITIAL_HEADER_TABLE(version_pre7_main_list); + +static Ecat_header_table_type +version_pre7_subhdr_table = INITIAL_HEADER_TABLE(version_pre7_subhdr_list); + +static Ecat_header_description_type version_pre7_description = +{ + &version_pre7_main_table, + &version_pre7_subhdr_table +}; + +static Ecat_header_description_type *ECAT_VER_PRE7 = &version_pre7_description; +
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/ecattominc.c @@ -0,0 +1,1446 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : ecattominc +@INPUT : argc, argv - command line arguments +@OUTPUT : (none) +@RETURNS : error status +@DESCRIPTION: Converts a CTI ECAT file to a minc format file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 3, 1996 (Peter Neelin) +@MODIFIED : $Log: ecattominc.c,v $ +@MODIFIED : Revision 1.1 1996-01-18 14:52:14 neelin +@MODIFIED : Initial revision +@MODIFIED : +@COPYRIGHT : + Copyright 1996 Peter Neelin, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ + +#ifndef lint +static char rcsid[]="$Header: /private-cvsroot/minc/conversion/ecattominc/ecattominc.c,v 1.1 1996-01-18 14:52:14 neelin Exp $"; +#endif + +#include <stdlib.h> +#include <stdio.h> +#include <string.h> +#include <math.h> +#include <limits.h> +#include <float.h> +#include <time.h> +#include <ParseArgv.h> +#include <time_stamp.h> +#include <minc.h> +#include <minc_def.h> +#include <ecat_file.h> + +/* Type declarations */ +typedef struct { + char *name; + char *values; +} ecat_header_data_type; + +typedef struct { + int nslices; + int low_slice; + int high_slice; + double scan_time; + double time_width; + double half_life; + double zstart; + double zstep; + char isotope[16]; + int image_xsize; + int image_ysize; + int ordered_frame; + int *ordered_slices; + char image_type[16]; +} frame_info_type; + +typedef struct { + int low_frame; + int high_frame; + int num_slices; + int max_nslices; + int max_xsize; + int max_ysize; + double zwidth; + double xstep; + double ystep; + double xwidth; + double ywidth; + int decay_corrected; + char img_units[16]; + char patient_name[40]; + char patient_sex[8]; + long patient_age; + char patient_birthdate[40]; + char study_id[40]; + char start_time[40]; + long start_year; + long start_month; + long start_day; + long start_hour; + long start_minute; + double start_seconds; + char tracer[40]; + char injection_time[40]; + long injection_hour; + long injection_minute; + double injection_seconds; + double injection_dose; + ecat_header_data_type *main_field_list; + int num_main_fields; + ecat_header_data_type *subhdr_field_list; + int num_subhdr_fields; +} general_info_type; + +typedef struct { + double sort_key; + void *sort_value; +} sort_type; + +/* Function declarations */ +void usage_error(char *progname); +int get_frame_info(Ecat_file *ecat_fp, int slice_range[2], + int num_frames, frame_info_type *frame_info, + general_info_type *general_info); +void sort_slices(int sort_over_time, int num_frames, + frame_info_type *frame_info, + general_info_type *general_info); +int sortcmp(const void *val1, const void *val2); +int setup_minc_file(int mincid, int write_byte_data, int copy_all_header, + int ndims, long count[], int num_frames, + frame_info_type *frame_info, + general_info_type *general_info, + char *blood_file); +int get_slice(Ecat_file *ecat_fp, int frame_num, int slice_num, + long *pixel_max, double *image_max, short *image, + frame_info_type *frame_info, + general_info_type *general_info); +int write_minc_slice(double scale, int write_byte_data, + int mincid, int icvid, + int ndims,long start[], long count[], + short *image, int image_xsize, int image_ysize, + long pixel_max, double image_max, + double scan_time, double time_width, double zpos); +double decay_correction(double scan_time, double measure_time, + double start_time, double half_life); +void CreateBloodStructures (int mincHandle, int bloodHandle); +void FillBloodStructures (int mincHandle, int bloodHandle); + +/* Constants */ +#define TRUE 1 +#define FALSE 0 +#define MAX_DIMS 4 +#define ECAT_ACTIVITY "ACTIVITY" +#define ECAT_CALIB_UNITS_CPS 3 +#define MM_PER_CM 10.0 +#define BECQUEREL_PER_MCURIE 3.7e7 +#define SECONDS_PER_HOUR 3600 +#define DEFAULT_RANGE INT_MIN + +/* Main program */ + +int main(int argc, char *argv[]) +{ + /* Variables for arguments */ + static int write_byte_data = TRUE; + static int clobber = FALSE; + static int verbose = TRUE; + static int decay_correct = TRUE; + static int slice_range[2] = {DEFAULT_RANGE, DEFAULT_RANGE}; + static int copy_all_header = FALSE; + static char *blood_file = NULL; + static int frame_range[2] = {DEFAULT_RANGE, DEFAULT_RANGE}; + + /* Argument option table */ + static ArgvInfo argTable[] = { + {"-byte", ARGV_CONSTANT, (char *) TRUE, (char *) &write_byte_data, + "Write out data as bytes (default)."}, + {"-short", ARGV_CONSTANT, (char *) FALSE, (char *) &write_byte_data, + "Write out data as short integers."}, + {"-decay_correct", ARGV_CONSTANT, (char *) TRUE, (char *) &decay_correct, + "Do decay correction on images (default)."}, + {"-nodecay_correct", ARGV_CONSTANT, (char *) FALSE, + (char *) &decay_correct, "Don't do decay correction."}, + {"-clobber", ARGV_CONSTANT, (char *) TRUE, (char *) &clobber, + "Overwrite existing file."}, + {"-noclobber", ARGV_CONSTANT, (char *) FALSE, (char *) &clobber, + "Don't overwrite existing file (default)."}, + {"-verbose", ARGV_CONSTANT, (char *) TRUE, (char *) &verbose, + "List files as they are converted (default)"}, + {"-quiet", ARGV_CONSTANT, (char *) FALSE, (char *) &verbose, + "Do not list files as they are converted."}, + {"-slices", ARGV_INT, (char *) 2, (char *) slice_range, + "Range of slices to copy (counting from 0)."}, + {"-frames", ARGV_INT, (char *) 2, (char *) frame_range, + "Range of frames to copy (counting from 0)."}, + {"-frame", ARGV_INT, (char *) 1, (char *) frame_range, + "Single frame to copy (counting from 0)."}, + {"-small_header", ARGV_CONSTANT, (char *) FALSE, + (char *) ©_all_header, + "Copy only basic header information (default)."}, + {"-all_header", ARGV_CONSTANT, (char *) TRUE, (char *) ©_all_header, + "Copy all header information."}, + {"-bloodfile", ARGV_STRING, (char *) 1, (char *) &blood_file, + "Insert blood data from this file."}, + {NULL, ARGV_END, NULL, NULL, NULL} + }; + + /* Other variables */ + char *pname; + char *mincfile; + char *ecat_filename; + int num_frames; + int num_bed_positions; + int sort_over_time; + frame_info_type *frame_info; + general_info_type *general_info; + long count[MAX_DIMS], start[MAX_DIMS]; + int ndims; + int mincid, icvid, varid; + int islice, iframe, i, slice_num, ifield, frame_num, low_frame, high_frame; + long pixel_max; + double image_max; + double scale; + short *image; + Ecat_file *ecat_fp; + int status; + char *tm_stamp; + double first_z, last_z, zstep; + + /* Get time stamp */ + tm_stamp = time_stamp(argc, argv); + + /* Check arguments */ + pname = argv[0]; + if (ParseArgv(&argc, argv, argTable, 0) || (argc != 3)) { + usage_error(pname); + } + + /* Check the slice range */ + if (slice_range[0] == DEFAULT_RANGE) { + slice_range[0] = 0; + slice_range[1] = INT_MAX; + } + else if (slice_range[1] == DEFAULT_RANGE) { + slice_range[1] = slice_range[0]; + } + if ((slice_range[0] < 0) || (slice_range[1] < 0) || + (slice_range[1] < slice_range[0])) { + (void) fprintf(stderr, "%s: Error in slice range: %d to %d.\n", + pname, slice_range[0], slice_range[1]); + exit(EXIT_FAILURE); + } + + /* Check the frame range */ + if (frame_range[0] == DEFAULT_RANGE) { + frame_range[0] = 0; + frame_range[1] = INT_MAX; + } + else if (frame_range[1] == DEFAULT_RANGE) { + frame_range[1] = frame_range[0]; + } + if ((frame_range[0] < 0) || (frame_range[1] < 0) || + (frame_range[1] < frame_range[0])) { + (void) fprintf(stderr, "%s: Error in frame range: %d to %d.\n", + pname, frame_range[0], frame_range[1]); + exit(EXIT_FAILURE); + } + + /* Get file names */ + ecat_filename = argv[1]; + mincfile = argv[2]; + + /* Open the ECAT file */ + if ((ecat_fp=ecat_open(ecat_filename))==NULL) { + (void) fprintf(stderr, "%s: Error opening file %s.\n", + pname, ecat_filename); + exit(EXIT_FAILURE); + } + + /* Get number of frames and bed positions to see if we are varying over + time or interleaving slices */ + num_frames = ecat_get_num_frames(ecat_fp); + num_bed_positions = ecat_get_num_bed_positions(ecat_fp); + sort_over_time = TRUE; + if ((num_frames > 1) && (num_bed_positions > 1)) { + (void) fprintf(stderr, + "%s: Cannot handle multiple frames and bed positions.\n", + pname); + exit(EXIT_FAILURE); + } + if (num_bed_positions > 1) { + num_frames = num_bed_positions; + sort_over_time = FALSE; + } + + /* Print log message */ + if (verbose) { + (void) fprintf(stderr, "Reading headers.\n"); + } + + /* Get frame range */ + low_frame = 0; + high_frame = num_frames - 1; + if (frame_range[0] > low_frame) + low_frame = frame_range[0]; + if (frame_range[1] < high_frame) + high_frame = frame_range[1]; + if (low_frame > high_frame) + low_frame = high_frame; + num_frames = high_frame - low_frame + 1; + + /* Read the files to get basic information */ + general_info=MALLOC(sizeof(*general_info)); + frame_info = MALLOC(num_frames * sizeof(*frame_info)); + general_info->low_frame = low_frame; + general_info->high_frame = high_frame; + status=get_frame_info(ecat_fp, slice_range, + num_frames, frame_info, general_info); + if (status >= 0) { + (void) fprintf(stderr, "%s: Error reading ECAT file %s on frame %d.\n", + pname, ecat_filename, status); + exit(EXIT_FAILURE); + } + + /* Allocate space for ordered slice list if sorting over z position */ + if (!sort_over_time) { + for (iframe=0; iframe<num_frames; iframe++) { + frame_info[iframe].ordered_slices = + MALLOC(frame_info[iframe].nslices * sizeof(int)); + } + } + + /* Sort the slices */ + sort_slices(sort_over_time, num_frames, + frame_info, general_info); + + /* Setup the minc file */ + if (sort_over_time && (num_frames>1)) { + ndims = 4; + count[0]=num_frames; + count[1]=general_info->max_nslices; + } + else { + ndims=3; + count[0]=general_info->num_slices; + } + count[ndims-1] = general_info->max_xsize; + count[ndims-2] = general_info->max_ysize; + mincid = micreate(mincfile, (clobber ? NC_CLOBBER : NC_NOCLOBBER)); + (void) miattputstr(mincid, NC_GLOBAL, MIhistory, tm_stamp); + icvid=setup_minc_file(mincid, write_byte_data, copy_all_header, + ndims, count, num_frames, + frame_info, general_info, blood_file); + if (icvid==MI_ERROR) { + (void) fprintf(stderr, + "%s: Error setting up minc file %s.\n", + pname, mincfile); + exit(EXIT_FAILURE); + } + + /* Initialize minc start and count vectors */ + for (i=0; i<ndims-2; i++) count[i]=1; + (void) miset_coords(ndims, (long) 0, start); + + /* Set up values for decay correction */ + scale = 1.0; + + /* Allocate an image buffer */ + image = MALLOC(general_info->max_xsize * general_info->max_ysize * + sizeof(short)); + + /* Print log message */ + if (verbose) { + (void) fprintf(stderr, "Copying frames:"); + (void) fflush(stderr); + } + + /* Loop through files */ + for (iframe=0; iframe<num_frames; iframe++) { + + /* Print log message */ + if (verbose) { + (void) fprintf(stderr, "."); + (void) fflush(stderr); + } + + /* Get decay correction (scan time is already measured from injection + time) */ + if (decay_correct && !general_info->decay_corrected && + (strcmp(frame_info[iframe].image_type, ECAT_ACTIVITY) == 0)) { + scale = decay_correction(frame_info[iframe].scan_time, + frame_info[iframe].time_width, + 0.0, + frame_info[iframe].half_life); + } + else { + scale = 1.0; + } + + /* Loop through slices */ + for (islice = 0; islice < frame_info[iframe].nslices; islice++) { + + /* Set the minc coordinates */ + if (sort_over_time && (num_frames>1)) { + start[0]=frame_info[iframe].ordered_frame; + start[1]=islice; + } + else if (sort_over_time) + start[0]=islice; + else + start[0]=frame_info[iframe].ordered_slices[islice]; + count[ndims-1] = frame_info[iframe].image_xsize; + count[ndims-2] = frame_info[iframe].image_ysize; + + /* Copy the slice */ + slice_num = islice + frame_info[iframe].low_slice; + frame_num = iframe + general_info->low_frame; + if (get_slice(ecat_fp, frame_num, slice_num, + &pixel_max, &image_max, image, + &frame_info[iframe], general_info) || + write_minc_slice(scale, write_byte_data, + mincid, icvid, ndims, start, count, image, + frame_info[iframe].image_xsize, + frame_info[iframe].image_ysize, + pixel_max, image_max, + frame_info[iframe].scan_time, + frame_info[iframe].time_width, + frame_info[iframe].zstep * + (double) slice_num + + frame_info[iframe].zstart)) { + (void) fprintf(stderr, + "%s: Error copying slice %d from frame %d.\n", + pname, slice_num, frame_num); + exit(EXIT_FAILURE); + } + + } /* End slice loop */ + + + } /* End frame loop */ + + /* Write out average z step and start for irregularly spaced slices */ + if ((ndims!=MAX_DIMS) && (num_frames>1)) { + start[0] = 0; + varid = ncvarid(mincid, MIzspace); + (void) mivarget1(mincid, varid, start, NC_DOUBLE, NULL, + &first_z); + start[0] = general_info->num_slices - 1; + (void) mivarget1(mincid, varid, start, NC_DOUBLE, NULL, + &last_z); + if (start[0] > 0) + zstep = (last_z - first_z) / ((double) start[0]); + else + zstep = 1.0; + (void) miattputdbl(mincid, varid, MIstep, zstep); + (void) miattputdbl(mincid, varid, MIstart, first_z); + } + + /* Close minc file */ + (void) miattputstr(mincid, ncvarid(mincid, MIimage), MIcomplete, MI_TRUE); + (void) miclose(mincid); + + /* Write out log message */ + if (verbose) { + (void) fprintf(stderr, "Done\n"); + (void) fflush(stderr); + } + + FREE(image); + if (!sort_over_time) { + for (iframe=0; iframe<num_frames; iframe++) { + FREE(frame_info[iframe].ordered_slices); + } + } + FREE(frame_info); + for (ifield=0; ifield < general_info->num_main_fields; ifield++) { + FREE(general_info->main_field_list[ifield].name); + FREE(general_info->main_field_list[ifield].values); + } + for (ifield=0; ifield < general_info->num_subhdr_fields; ifield++) { + FREE(general_info->subhdr_field_list[ifield].name); + FREE(general_info->subhdr_field_list[ifield].values); + } + FREE(general_info->main_field_list); + FREE(general_info->subhdr_field_list); + FREE(general_info); + + exit(EXIT_SUCCESS); + +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : usage_error +@INPUT : progname - program name +@OUTPUT : (none) +@RETURNS : (nothing) +@DESCRIPTION: Prints a usage error message and exits. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 3, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +void usage_error(char *progname) +{ + (void) fprintf(stderr, + "\nUsage: %s [<options>] <infile> <outfile.mnc>\n", progname); + (void) fprintf(stderr, + " %s [-help]\n\n", progname); + + exit(EXIT_FAILURE); +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : get_frame_info +@INPUT : ecat_fp - file pointer for ecat file + slice_range - 2-component array giving range of slices + num_frames - number of frames +@OUTPUT : frame_info - array of structures containing information + about each frame. + general_info - general information about the file. +@RETURNS : (-1) if no error occurs, otherwise, the index + of the first frame that could not be read. +@DESCRIPTION: Reads information for frame in the ECAT file +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 4, 1996 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +int get_frame_info(Ecat_file *ecat_fp, int slice_range[2], + int num_frames, frame_info_type *frame_info, + general_info_type *general_info) +{ + static char *the_months[]= + {NULL, "Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", + "Aug", "Sep", "Oct", "Nov", "Dec"}; + int *num_slices, *max_nslices, *max_xsize, *max_ysize; + int start_day, start_month, start_year; + int start_hour, start_minute, start_seconds; + frame_info_type *fip; + int lvalue; + double fvalue; + char svalue[ECAT_MAX_STRING_LENGTH]; + int iframe, ifield, imult, num_fields, iheader, curframe; + int length, newlength, multiplicity; + struct tm *tm_ptr; + ecat_header_data_type *field_list; + Ecat_field_name field; + char *description; + time_t the_time; + + /* Initialize number of slices */ + num_slices = &(general_info->num_slices); + max_nslices = &(general_info->max_nslices); + max_xsize = &(general_info->max_xsize); + max_ysize = &(general_info->max_ysize); + *num_slices = 0; + *max_nslices = 0; + *max_xsize = *max_ysize = 0; + general_info->decay_corrected = FALSE; + + /* Loop through files, reading information */ + for (iframe=0; iframe<num_frames; iframe++) { + + /* Get pointer to frame_info structure */ + fip = &frame_info[iframe]; + curframe = iframe + general_info->low_frame; + + /* Get number of slices */ + fip->low_slice = 0; + fip->high_slice = ecat_get_num_planes(ecat_fp) - 1; + if (slice_range[0] > fip->low_slice) + fip->low_slice = slice_range[0]; + if (slice_range[1] < fip->high_slice) + fip->high_slice = slice_range[1]; + if (fip->low_slice > fip->high_slice) + fip->low_slice = fip->high_slice; + fip->nslices = fip->high_slice - fip->low_slice + 1; + *num_slices += fip->nslices; + if (fip->nslices > *max_nslices) + *max_nslices = fip->nslices; + + /* Get image width */ + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_X_Dimension, 0, + &lvalue, NULL, NULL)) return curframe; + fip->image_xsize = lvalue; + if (lvalue > *max_xsize) *max_ysize = lvalue; + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Y_Dimension, 0, + &lvalue, NULL, NULL)) return curframe; + fip->image_ysize = lvalue; + if (lvalue > *max_xsize) *max_xsize = lvalue; + + /* Get frame start time (in seconds) */ + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Frame_Start_Time, 0, + &lvalue, NULL, NULL)) return curframe; + fip->scan_time = (double) lvalue / 1000.0; + + /* Get length of frame (in seconds) */ + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Frame_Duration, 0, + &lvalue, NULL, NULL)) return curframe; + fip->time_width = (double) lvalue / 1000.0; + + /* Get scan type */ + if (ecat_get_main_value(ecat_fp, ECAT_Calibration_Units, 0, + &lvalue, NULL, NULL)) return curframe; + if (lvalue == ECAT_CALIB_UNITS_CPS) { + (void) strcpy(fip->image_type, ECAT_ACTIVITY); + } + else { + (void) strcpy(fip->image_type, ""); + } + + /* Get isotope and half-life */ + if (ecat_get_main_value(ecat_fp, ECAT_Isotope_Name, 0, + NULL, NULL, fip->isotope)) + return curframe; + if (ecat_get_main_value(ecat_fp, ECAT_Isotope_Halflife, 0, + NULL, &fip->half_life, NULL)) + return curframe; + + /* Get z start and step (correct start for non-zero first slice */ + if (ecat_get_main_value(ecat_fp, ECAT_Plane_Separation, 0, + NULL, &fip->zstep, NULL)) return curframe; + fip->zstep *= -MM_PER_CM; + if (ecat_get_num_bed_positions(ecat_fp) <= 1) { + if (ecat_get_main_value(ecat_fp, ECAT_Init_Bed_Position, 0, + NULL, &fip->zstart, NULL)) return curframe; + } + else { + if (ecat_get_main_value(ecat_fp, ECAT_Bed_Position, curframe, + NULL, &fip->zstart, NULL)) return curframe; + } + fip->zstart *= -MM_PER_CM; + fip->zstart += fip->low_slice * fip->zstep; + + /* Check to see if file has been decay corrected */ + if (!ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Decay_Corr_Fctr, 0, + NULL, &fvalue, NULL) && fvalue != 1.0) { + general_info->decay_corrected = TRUE; + } + + /* Get general information from first frame */ + if (iframe==0) { + + /* Get pixel sizes */ + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_X_Pixel_Size, 0, + NULL, &general_info->xstep, NULL)) + return curframe; + if (ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Y_Pixel_Size, 0, + NULL, &general_info->ystep, NULL)) + return curframe; + general_info->xstep *= MM_PER_CM; + general_info->ystep *= MM_PER_CM; + + /* Get resolution in each direction (or zero if not found) */ + general_info->xwidth = general_info->ywidth = + general_info->zwidth = -1.0; + (void) ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_X_Resolution, + 0, NULL, &general_info->xwidth, NULL); + (void) ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Y_Resolution, + 0, NULL, &general_info->ywidth, NULL); + (void) ecat_get_subhdr_value(ecat_fp, curframe, 0, ECAT_Z_Resolution, + 0, NULL, &general_info->zwidth, NULL); + general_info->xwidth *= MM_PER_CM; + general_info->ywidth *= MM_PER_CM; + general_info->zwidth *= MM_PER_CM; + + /* Get image range and units */ + if (ecat_get_main_value(ecat_fp, ECAT_Calibration_Units, 0, + &lvalue, NULL, NULL)) return curframe; + if (lvalue == ECAT_CALIB_UNITS_CPS) { + (void) strcpy(general_info->img_units, "nCi/cc"); + } + else { + (void) strcpy(general_info->img_units, ""); + } + + /* Get patient information */ + if (ecat_get_main_value(ecat_fp, ECAT_Patient_Name, 0, + NULL, NULL, general_info->patient_name)) + return curframe; + if (ecat_get_main_value(ecat_fp, ECAT_Patient_Sex, 0, + NULL, NULL, general_info->patient_sex)) + return curframe; + switch (general_info->patient_sex[0]) { + case 1: + case 'M': + (void) strcpy(general_info->patient_sex, MI_MALE); + break; + case 2: + case 'F': + (void) strcpy(general_info->patient_sex, MI_FEMALE); + break; + default: + (void) strcpy(general_info->patient_sex, MI_OTHER); + break; + } + if (ecat_get_main_value(ecat_fp, ECAT_Patient_Age, 0, + &lvalue, NULL, NULL)) + general_info->patient_age = -1; + else + general_info->patient_age = lvalue; + if (!ecat_get_main_value(ecat_fp, ECAT_Patient_Birth_Date, 0, + &lvalue, NULL, NULL)) { + /* Try to get the right birthday by adding half a day, using + UTC and rounding down. This works because field stores + birthday at 0:00 converted using the local scanner + timezone. */ + the_time = (time_t) (lvalue + 12 * SECONDS_PER_HOUR); + tm_ptr = gmtime(&the_time); + (void) sprintf(general_info->patient_birthdate, + "%d-%s-%d", tm_ptr->tm_mday, + the_months[tm_ptr->tm_mon+1], + tm_ptr->tm_year+1900); + } + else { + general_info->patient_birthdate[0] = '\0'; + } + + /* Get study information */ + if (ecat_get_main_value(ecat_fp, ECAT_Study_Type, 0, + NULL, NULL, general_info->study_id)) + return curframe; + if (!ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Time, 0, + &lvalue, NULL, NULL)) { + the_time = (time_t) lvalue; + tm_ptr = localtime(&the_time); + general_info->start_day = tm_ptr->tm_mday; + general_info->start_month = tm_ptr->tm_mon + 1; + general_info->start_year = tm_ptr->tm_year; + general_info->start_hour = tm_ptr->tm_hour; + general_info->start_minute = tm_ptr->tm_min; + general_info->start_seconds = tm_ptr->tm_sec; + } + else { + if (ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Day, 0, + &start_day, NULL, NULL) || + ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Month, 0, + &start_month, NULL, NULL) || + ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Year, 0, + &start_year, NULL, NULL) || + ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Hour, 0, + &start_hour, NULL, NULL) || + ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Minute, 0, + &start_minute, NULL, NULL) || + ecat_get_main_value(ecat_fp, ECAT_Scan_Start_Second, 0, + &start_seconds, NULL, NULL)) { + return curframe; + } + general_info->start_day = start_day; + general_info->start_month = start_month; + general_info->start_year = start_year; + general_info->start_hour = start_hour; + general_info->start_minute = start_minute; + general_info->start_seconds = start_seconds; + } + (void) sprintf(general_info->start_time, "%d-%s-%d %d:%d:%d", + (int) general_info->start_day, + the_months[general_info->start_month], + (int) general_info->start_year+1900, + (int) general_info->start_hour, + (int) general_info->start_minute, + (int) general_info->start_seconds); + if (ecat_get_main_value(ecat_fp, ECAT_Radiopharmaceutical, 0, + NULL, NULL, general_info->tracer)) + return curframe; + if (!ecat_get_main_value(ecat_fp, ECAT_Dose_Start_Time, 0, + &lvalue, NULL, NULL) && (lvalue != 0)) { + the_time = (time_t) lvalue; + (void) strcpy(general_info->injection_time, + ctime(&the_time)); + tm_ptr = localtime(&the_time); + general_info->injection_hour = tm_ptr->tm_hour; + general_info->injection_minute = tm_ptr->tm_min; + general_info->injection_seconds = tm_ptr->tm_sec; + } + else { + general_info->injection_time[0] = '\0'; + } + if (!ecat_get_main_value(ecat_fp, ECAT_Dosage, 0, + NULL, &general_info->injection_dose, NULL)) { + general_info->injection_dose /= BECQUEREL_PER_MCURIE; + } + else { + general_info->injection_dose = -1.0; + } + + /* Get list of header values */ + for (iheader=0; iheader < 2; iheader++) { + + /* Get space for first field */ + field_list = MALLOC(sizeof(*field_list)); + + /* Loop through fields */ + ifield = 0; + num_fields = 0; + do { + + /* Get next field */ + if (iheader == 0) { + field = ecat_list_main(ecat_fp, ifield); + description = + ecat_get_main_field_description(ecat_fp, field); + multiplicity = + ecat_get_main_field_length(ecat_fp, field); + } + else { + field = ecat_list_subhdr(ecat_fp, ifield); + description = + ecat_get_subhdr_field_description(ecat_fp, field); + multiplicity = + ecat_get_subhdr_field_length(ecat_fp, field); + } + ifield++; + + /* Check for end of list */ + if ((field == ECAT_No_Field) || (description == NULL) || + (multiplicity <= 0)) + continue; + + /* Save the description */ + field_list[num_fields].name = strdup(description); + + /* Get space for the values */ + field_list[num_fields].values = NULL; + length = 0; + + /* Loop through multiplicity */ + for (imult=0; imult < multiplicity; imult++) { + + /* Get value */ + svalue[0] = '\0'; + if (iheader == 0) { + (void) ecat_get_main_value(ecat_fp, field, imult, + NULL, NULL, svalue); + } + else { + (void) ecat_get_subhdr_value(ecat_fp, curframe, 0, + field, imult, + NULL, NULL, svalue); + } + + /* Save it */ + if (field_list[num_fields].values == NULL) { + field_list[num_fields].values = strdup(svalue); + length = strlen(svalue); + } + else { + newlength = length + strlen(svalue) + 1; + field_list[num_fields].values = + REALLOC(field_list[num_fields].values, + (size_t) newlength + 1); + field_list[num_fields].values[length] = '\n'; + (void) strcpy(&field_list[num_fields].values[length+1], + svalue); + length = newlength; + } + + } /* End of loop over multiplicity */ + + + /* Get space for next field */ + field_list = REALLOC(field_list, + (num_fields+2) * sizeof(field_list[0])); + + /* Increment counter */ + num_fields++; + + } while (field != ECAT_No_Field); + + /* Save the list */ + if (iheader == 0) { + general_info->main_field_list=field_list; + general_info->num_main_fields=num_fields; + } + else { + general_info->subhdr_field_list=field_list; + general_info->num_subhdr_fields=num_fields; + } + + } /* Loop over headers */ + + } /* If first file */ + + } /* Loop over files */ + + return -1; +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : sort_slices +@INPUT : sort_over_time - boolean indicating whether sort should be + over time or z position. + num_frames - number of frames + frame_info - array of frame information. + general_info - general information about files. +@OUTPUT : frame_info - modified to give slice ordering information. +@RETURNS : (nothing) +@DESCRIPTION: Sorts the slices for the output file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 12, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +void sort_slices(int sort_over_time, int num_frames, + frame_info_type *frame_info, + general_info_type *general_info) +{ + int iframe, islice, isort, num_sort, slice_num; + /* Variables for sorting */ + sort_type *sort_array; + struct { + int file; + int slice; + } *slice_ptr; + frame_info_type *file_ptr; + + /* Allocate array for sorting */ + num_sort = (sort_over_time ? num_frames : general_info->num_slices); + sort_array = MALLOC (num_sort * sizeof(*sort_array)); + + /* Are we sorting over time or z position */ + if (sort_over_time) { + /* Go through the files */ + for (iframe=0; iframe<num_frames; iframe++) { + sort_array[iframe].sort_key = frame_info[iframe].scan_time; + sort_array[iframe].sort_value = &frame_info[iframe]; + } + } + /* Otherwise we are sorting over z position */ + else { + /* Go through the files and slices */ + isort=0; + for (iframe=0; iframe<num_frames; iframe++) { + for (islice = 0; islice < frame_info[iframe].nslices; islice++) { + slice_num = islice + frame_info[iframe].low_slice; + sort_array[isort].sort_key = frame_info[iframe].zstart + + slice_num * frame_info[iframe].zstep; + slice_ptr = MALLOC(sizeof(*slice_ptr)); + slice_ptr->file = iframe; + slice_ptr->slice = islice; + sort_array[isort].sort_value = slice_ptr; + isort++; + } + } + } + + /* Sort the slices */ + qsort(sort_array, num_sort, sizeof(*sort_array), sortcmp); + + /* Loop through sorted list */ + for (isort=0; isort<num_sort; isort++) { + if (sort_over_time) { + file_ptr=sort_array[isort].sort_value; + file_ptr->ordered_frame = isort; + } + else { + slice_ptr = sort_array[isort].sort_value; + iframe = slice_ptr->file; + islice = slice_ptr->slice; + frame_info[iframe].ordered_slices[islice] = isort; + FREE(slice_ptr); + } + } + + /* Free the sorting array */ + FREE(sort_array); + + return; +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : sortcmp +@INPUT : val1 - first value + val2 - second value +@OUTPUT : (none) +@RETURNS : 0 if values are the same, -1 if val1->sort_key < val2->sort_key + and +1 if val1->sort_key > val2->sort_key. +@DESCRIPTION: Compares two double precision values. If they are the same, + then return 0. If val1 < val2, return -1. If val1 > val2, + return +1. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 12, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +int sortcmp(const void *val1, const void *val2) +{ + + if (((sort_type *)val1)->sort_key < + ((sort_type *)val2)->sort_key) + return -1; + else if (((sort_type *)val1)->sort_key > + ((sort_type *)val2)->sort_key) + return 1; + else return 0; + +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : setup_minc_file +@INPUT : mincid - id of minc file + write_byte_data - boolean indicating whether data should be + written as bytes (TRUE) or shorts (FALSE). + copy_all_header - boolean indicating whether all of the + header information should be copied or not. + ndims - number of dimensions for minc file + count - lengths of dimensions minc file + num_frames - number of frames. + frame_info - array of information about frames. + general_info - general information about file. + blood_file - name of blood file containing data to include. +@OUTPUT : (nothing) +@RETURNS : Image conversion variable id or MI_ERROR if an error occurs. +@DESCRIPTION: Initializes the header of the minc file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 12, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +int setup_minc_file(int mincid, int write_byte_data, int copy_all_header, + int ndims, long count[], int num_frames, + frame_info_type *frame_info, + general_info_type *general_info, + char *blood_file) +{ + static char *dim_names_array[]={MItime, MIzspace, MIyspace, MIxspace}; + char **dim_names; + static char *dimwidth_names_array[]={ + MItime_width, MIzspace_width, MIyspace_width, MIxspace_width}; + char **dimwidth_names; + int dim[MAX_DIMS]; + int img, imgmax, imgmin, dimvarid, widvarid, icv, varid, ecat_var; + int bloodid; + int idim, ifield, num_fields, iheader; + double vrange[2]; + char varname[MAX_NC_NAME]; + ecat_header_data_type *field_list; + double dimwidths[MAX_DIMS]; + + /* Set up dimension arrays for looping */ + dim_names = dim_names_array + MAX_DIMS - ndims; + dimwidth_names = dimwidth_names_array + MAX_DIMS - ndims; + for (idim=0; idim < ndims; idim++) { + switch (idim + MAX_DIMS - ndims) { + case 0: dimwidths[idim] = 1.0; break; + case 1: dimwidths[idim] = general_info->zwidth; break; + case 2: dimwidths[idim] = general_info->ywidth; break; + case 3: dimwidths[idim] = general_info->xwidth; break; + } + } + + /* Create the dimensions */ + for (idim=0; idim<ndims; idim++) { + dim[idim]=ncdimdef(mincid, dim_names[idim], count[idim]); + dimvarid=micreate_std_variable(mincid, dim_names[idim], NC_DOUBLE, + ((idim==0) && (num_frames>1)) ? 1 : 0, + &dim[idim]); + if (dimwidths[idim] > 0) { + widvarid=micreate_std_variable(mincid, dimwidth_names[idim], + NC_DOUBLE, + ((strcmp(dim_names[idim], MItime)==0) + ? 1 : 0), &dim[idim]); + } + else { + widvarid = MI_ERROR; + } + + /* Add attributes to the dimension variables */ + if (strcmp(dim_names[idim], MIzspace)==0) { + /* Write out step and start. We will rewrite this for irregularly + spaced files */ + (void) miattputdbl(mincid, dimvarid, MIstep, + frame_info[0].zstep); + (void) miattputdbl(mincid, dimvarid, MIstart, + frame_info[0].zstart); + (void) miattputstr(mincid, dimvarid, MIunits, "mm"); + (void) miattputstr(mincid, dimvarid, MIspacetype, MI_NATIVE); + if (widvarid != MI_ERROR) { + (void) miattputdbl(mincid, widvarid, MIwidth, + (double) general_info->zwidth); + (void) miattputstr(mincid, widvarid, MIunits, "mm"); + (void) miattputstr(mincid, widvarid, MIfiltertype, MI_GAUSSIAN); + } + } + else if (strcmp(dim_names[idim], MIyspace)==0) { + (void) miattputstr(mincid, dimvarid, MIunits, "mm"); + (void) miattputdbl(mincid, dimvarid, MIstep, + (double) general_info->ystep); + (void) miattputstr(mincid, dimvarid, MIspacetype, MI_NATIVE); + if (widvarid != MI_ERROR) { + (void) miattputdbl(mincid, widvarid, MIwidth, + (double) general_info->ywidth); + (void) miattputstr(mincid, widvarid, MIfiltertype, MI_GAUSSIAN); + } + } + else if (strcmp(dim_names[idim], MIxspace)==0) { + (void) miattputstr(mincid, dimvarid, MIunits, "mm"); + (void) miattputdbl(mincid, dimvarid, MIstep, + (double) general_info->xstep); + (void) miattputstr(mincid, dimvarid, MIspacetype, MI_NATIVE); + if (widvarid != MI_ERROR) { + (void) miattputdbl(mincid, widvarid, MIwidth, + (double) general_info->xwidth); + (void) miattputstr(mincid, widvarid, MIfiltertype, MI_GAUSSIAN); + } + } + else if (strcmp(dim_names[idim], MItime)==0) { + (void) miattputstr(mincid, dimvarid, MIunits, "seconds"); + (void) miattputstr(mincid, widvarid, MIunits, "seconds"); + } + } + + /* Create the image variable */ + if (write_byte_data) { + img=micreate_std_variable(mincid, MIimage, NC_BYTE, ndims, dim); + (void) miattputstr(mincid, img, MIsigntype, MI_UNSIGNED); + vrange[0]=0; vrange[1]=255; + } + else { + img=micreate_std_variable(mincid, MIimage, NC_SHORT, ndims, dim); + (void) miattputstr(mincid, img, MIsigntype, MI_SIGNED); + vrange[0] = -32000; + vrange[1] = 32000; + } + (void) ncattput(mincid, img, MIvalid_range, NC_DOUBLE, 2, vrange); + (void) miattputstr(mincid, img, MIcomplete, MI_FALSE); + + /* Create the image max and min variables */ + imgmax=micreate_std_variable(mincid, MIimagemax, NC_DOUBLE, ndims-2, dim); + imgmin=micreate_std_variable(mincid, MIimagemin, NC_DOUBLE, ndims-2, dim); + (void) miattputstr(mincid, imgmax, MIunits, + general_info->img_units); + (void) miattputstr(mincid, imgmin, MIunits, + general_info->img_units); + + /* Create the image conversion variable */ + icv=miicv_create(); + (void) miicv_setint(icv, MI_ICV_TYPE, NC_SHORT); + + /* Save patient info */ + varid = micreate_group_variable(mincid, MIpatient); + (void) miattputstr(mincid, varid, MIfull_name, + general_info->patient_name); + (void) miattputstr(mincid, varid, MIsex, + general_info->patient_sex); + if (general_info->patient_age > 0) + (void) ncattput(mincid, varid, MIage, NC_LONG, 1, + &general_info->patient_age); + if ((int) strlen(general_info->patient_birthdate) > 0) + (void) miattputstr(mincid, varid, MIbirthdate, + general_info->patient_birthdate); + + /* Save study info */ + varid = micreate_group_variable(mincid, MIstudy); + (void) miattputstr(mincid, varid, MImodality, MI_PET); + (void) miattputstr(mincid, varid, MImanufacturer, "CTI"); + (void) miattputstr(mincid, varid, MIstudy_id, + general_info->study_id); + (void) miattputstr(mincid, varid, MIstart_time, + general_info->start_time); + (void) ncattput(mincid, varid, MIstart_year, NC_LONG, 1, + &general_info->start_year); + (void) ncattput(mincid, varid, MIstart_month, NC_LONG, 1, + &general_info->start_month); + (void) ncattput(mincid, varid, MIstart_day, NC_LONG, 1, + &general_info->start_day); + (void) ncattput(mincid, varid, MIstart_hour, NC_LONG, 1, + &general_info->start_hour); + (void) ncattput(mincid, varid, MIstart_minute, NC_LONG, 1, + &general_info->start_minute); + (void) ncattput(mincid, varid, MIstart_seconds, NC_DOUBLE, 1, + &general_info->start_seconds); + + /* Save acquisition info */ + varid = micreate_group_variable(mincid, MIacquisition); + (void) miattputstr(mincid, varid, MIradionuclide, + frame_info[0].isotope); + if (frame_info[0].half_life > 0.0) { + (void) miattputdbl(mincid, varid, MIradionuclide_halflife, + (double) frame_info[0].half_life); + } + (void) miattputstr(mincid, varid, MItracer, + general_info->tracer); + if ((int) strlen(general_info->injection_time) > 0) { + (void) miattputstr(mincid, varid, MIinjection_time, + general_info->injection_time); + (void) ncattput(mincid, varid, MIinjection_hour, NC_LONG, 1, + &general_info->injection_hour); + (void) ncattput(mincid, varid, MIinjection_minute, NC_LONG, 1, + &general_info->injection_minute); + (void) ncattput(mincid, varid, MIinjection_seconds, NC_DOUBLE, 1, + &general_info->injection_seconds); + } + if (general_info->injection_dose > 0.0) { + (void) ncattput(mincid, varid, MIinjection_dose, NC_DOUBLE, 1, + &general_info->injection_dose); + (void) miattputstr(mincid, varid, MIdose_units, "mCurie"); + } + + /* If we want all of the values from the header, get them */ + if (copy_all_header) { + + for (iheader=0; iheader < 2; iheader++) { + + /* Set up values for either header */ + if (iheader == 0) { + (void) strcpy(varname, "ecat-main"); + field_list = general_info->main_field_list; + num_fields = general_info->num_main_fields; + } + else { + (void) strcpy(varname, "ecat-subhdr"); + field_list = general_info->subhdr_field_list; + num_fields = general_info->num_subhdr_fields; + } + + /* Create a variable for ECAT fields */ + ecat_var = ncvardef(mincid, varname, NC_LONG, 0, NULL); + (void) miattputstr(mincid, ecat_var, MIvartype, MI_GROUP); + (void) miattputstr(mincid, ecat_var, MIvarid, "MNI ECAT variable"); + (void) miadd_child(mincid, ncvarid(mincid, MIrootvariable), ecat_var); + + /* Loop through fields */ + for (ifield=0; ifield < num_fields; ifield++){ + (void) miattputstr(mincid, ecat_var, field_list[ifield].name, + field_list[ifield].values); + } + + } /* Loop over headers */ + + } /* If copy_all_header */ + + /* Open the blood file and create the variables if needed */ + if (blood_file != NULL) { + bloodid = ncopen(blood_file, NC_NOWRITE); + CreateBloodStructures(mincid, bloodid); + } + + /* Attach the icv */ + (void) ncsetfill(mincid, NC_NOFILL); + (void) ncendef(mincid); + (void) miicv_attach(icv, mincid, img); + + /* Copy the blood data */ + if (blood_file != NULL) { + FillBloodStructures(mincid, bloodid); + ncclose(bloodid); + } + + return icv; +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : get_slice +@INPUT : ecat_fp - file pointer for file + frame_num - number of frame + slice_num - slice to copy + frame_info - information on frame + general_info - general file information +@OUTPUT : pixel_max - maximum pixel value + image_max - real value to which pixel_max corresponds + image - the image +@RETURNS : Returns TRUE if an error occurs. +@DESCRIPTION: Gets an image from the file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 20, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +int get_slice(Ecat_file *ecat_fp, int frame_num, int slice_num, + long *pixel_max, double *image_max, short *image, + frame_info_type *frame_info, + general_info_type *general_info) + /* ARGSUSED */ +{ + long npix, ix, iy, y_offset, off1, off2; + int pmax; + short temp; + double scale, global_scale; + + /* Get the image from the file */ + if (ecat_get_image(ecat_fp, frame_num, slice_num, image)) return TRUE; + + /* Flip the image to give positive x & y axes */ + npix = frame_info->image_xsize * frame_info->image_ysize; + for (iy=0; iy<frame_info->image_ysize/2; iy++) { + y_offset = iy * frame_info->image_xsize; + for (ix=0; ix<frame_info->image_xsize; ix++) { + off1 = y_offset + ix; + off2 = npix - off1 - 1; + temp = image[off1]; + image[off1] = image[off2]; + image[off2] = temp; + } + } + + /* Get image and pixel max */ + if (ecat_get_subhdr_value(ecat_fp, frame_num, slice_num, + ECAT_Image_Max, 0, &pmax, NULL, NULL) || + ecat_get_subhdr_value(ecat_fp, frame_num, slice_num, + ECAT_Scale_Factor, 0, NULL, &scale, NULL)) + return TRUE; + global_scale = 1.0; + (void) ecat_get_main_value(ecat_fp, ECAT_Calibration_Factor, 0, + NULL, &global_scale, NULL); + if (global_scale <= 0.0) global_scale = 1.0; + *pixel_max = pmax; + *image_max = (double) pmax * scale * global_scale; + + return FALSE; +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : write_minc_slice +@INPUT : scale - scale for decay correcting image + write_byte_data - boolean indicating whether data should be + written as bytes (TRUE) or shorts (FALSE). + mincid - id of minc file + icvid - id of image conversion variable + start - coordinate of slice in minc file + count - edge lengths of image to write in minc file + image - pointer to image buffer + image_xsize, image_ysize - dimensions of image + pixel_max - maximum pixel value + image_max - real value to which pixel_max corresponds + scan_time - time of slice + time_width - time width of slice + zpos - z position of slice + frame_info - information on frame + general_info - general file information +@OUTPUT : (nothing) +@RETURNS : Returns TRUE if an error occurs. +@DESCRIPTION: Writes out the image to the minc file. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 12, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +int write_minc_slice(double scale, int write_byte_data, + int mincid, int icvid, + int ndims,long start[], long count[], + short *image, int image_xsize, int image_ysize, + long pixel_max, double image_max, + double scan_time, double time_width, double zpos) + /*ARGSUSED*/ +{ + double pixmin, pixmax; + long ipix, npix; + double maximum, minimum; + + /* Search for pixel max and min */ + npix = image_xsize * image_ysize; + pixmin = pixmax = image[0]; + for (ipix=1; ipix<npix; ipix++) { + if (image[ipix]>pixmax) pixmax = image[ipix]; + if (image[ipix]<pixmin) pixmin = image[ipix]; + } + + /* Get image max and min */ + maximum = image_max * pixmax / (double) pixel_max; + minimum = image_max * pixmin / (double) pixel_max; + + /* Change valid range on icv */ + (void) miicv_detach(icvid); + (void) miicv_setdbl(icvid, MI_ICV_VALID_MAX, (double) pixmax); + (void) miicv_setdbl(icvid, MI_ICV_VALID_MIN, (double) pixmin); + (void) miicv_attach(icvid, mincid, ncvarid(mincid, MIimage)); + + /* Calculate real max and min */ + maximum = maximum * scale; + minimum = minimum * scale; + + /* Write out the image */ + (void) miicv_put(icvid, start, count, image); + + /* Write out image max and min */ + (void) ncvarput1(mincid, ncvarid(mincid, MIimagemax), start, &maximum); + (void) ncvarput1(mincid, ncvarid(mincid, MIimagemin), start, &minimum); + + /* Write out time and z position */ + if (ndims==MAX_DIMS) { + (void) mivarput1(mincid, ncvarid(mincid, MItime), start, + NC_DOUBLE, NULL, &scan_time); + (void) mivarput1(mincid, ncvarid(mincid, MItime_width), start, + NC_DOUBLE, NULL, &time_width); + } + (void) mivarput1(mincid, ncvarid(mincid, MIzspace), &start[ndims-3], + NC_DOUBLE, NULL, &zpos); + + return FALSE; +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : decay_correction +@INPUT : scan_time - time of beginning of sample + measure_time - length of sample + start_time - time to which we should decay correct + half_life - half life of isotope +@OUTPUT : (nothing) +@RETURNS : decay correction scaling factor +@DESCRIPTION: Calculates the decay correction needed to get equivalent counts + at start_time. Correction assumes constant activity (without + decay) over measure_time. +@METHOD : +@GLOBALS : +@CALLS : +@CREATED : January 21, 1993 (Peter Neelin) +@MODIFIED : +---------------------------------------------------------------------------- */ +double decay_correction(double scan_time, double measure_time, + double start_time, double half_life) +{ + double mean_life; + double measure_correction; + double decay; + + /* Check for negative half_life and calculate mean life */ + if (half_life <= 0.0) return 1.0; + mean_life = half_life/ log(2.0); + + /* Normalize scan time and measure_time */ + scan_time = (scan_time - start_time) / mean_life; + measure_time /= mean_life; + + /* Calculate correction for decay over measuring time (assuming a + constant activity). Check for possible rounding errors. */ + if ((measure_time*measure_time/2.0) < DBL_EPSILON) { + measure_correction = 1.0 - measure_time/2.0; + } + else { + measure_correction = (1.0 - exp(-measure_time)) / fabs(measure_time); + } + + /* Calculate decay */ + decay = exp(-scan_time) * measure_correction; + if (decay<=0.0) decay = DBL_MAX; + else decay = 1.0/decay; + + return decay; +} +
new file mode 100644 --- /dev/null +++ b/conversion/ecattominc/insertblood.c @@ -0,0 +1,70 @@ +#include <stdio.h> +#include <stdlib.h> +#include <minc.h> + +#define MIbloodroot "blood_analysis" + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : CreateBloodStructures +@INPUT : mincHandle -> a handle for an open MINC file. This file should + be open for writing, but not in redefinition + mode. + bloodHandle -> a handle for an open BNC file. This file should + be open for reading. +@OUTPUT : none +@RETURNS : void +@DESCRIPTION: Copies all variable definitions (with attributes) from the BNC + file to the MINC file. The appropriate dimensions are also + copied. +@METHOD : none. Just muddled through. +@GLOBALS : none +@CALLS : micopy_all_var_defs (MINC library) + miadd_child (MINC library) +@CREATED : May 30, 1994 by MW +@MODIFIED : +---------------------------------------------------------------------------- */ +void CreateBloodStructures (int mincHandle, int bloodHandle) +{ + int mincRoot; + int bloodRoot; + + /* + * Copy all the variables with their attributes. + */ + + (void) micopy_all_var_defs (bloodHandle, mincHandle, 0, NULL); + + /* + * Make the blood analysis root variable a child of + * the MINC root variable. + */ + + mincRoot = ncvarid (mincHandle, MIrootvariable); + bloodRoot = ncvarid (mincHandle, MIbloodroot); + (void) miadd_child (mincHandle, mincRoot, bloodRoot); +} + + +/* ----------------------------- MNI Header ----------------------------------- +@NAME : FillBloodStructures +@INPUT : mincHandle -> a handle for an open MINC file. This file should + be open for writing, but not in redefinition + mode. + bloodHandle -> a handle for an open BNC file. This file should + be open for reading. +@OUTPUT : none +@RETURNS : void +@DESCRIPTION: Copies all variable values from the BNC file to the MINC file. + The variable themselves should already exist in the MINC file + (see CreateBloodStructures). +@METHOD : none. Just muddled through. +@GLOBALS : none +@CALLS : micopy_all_var_values (MINC library) +@CREATED : May 30, 1994 by MW +@MODIFIED : +---------------------------------------------------------------------------- */ +void FillBloodStructures (int mincHandle, int bloodHandle) +{ + (void) micopy_all_var_values (bloodHandle, mincHandle, 0, NULL); +}