/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * filename: m-conc.c * * * * UTIL C-source: Medical Image Conversion Utility * * * * purpose : Read and Write Concorde microPET format files * * * * project : (X)MedCon by Erik Nolf * * * * Functions : MdcLoadPlaneCONC() - Load in 1 plane of file * * MdcLoadHeaderCONC() - Load in Concorde header info * * MdcLoadCONC() - Load Concorde file * * MdcSavePlaneCONC - Writeout one plane of file * * MdcSaveHeaderCONC - Writeout Concorde header info * * MdcSaveCONC() - Save Concorde file * * MdcCheckCONC() - Check for Concorde format * * MdcReadCONC() - Read Concorde file * * MdcWriteCONC() - Write Concorde file * * * * * * Author : Andy Loening * * * * Credits : * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ /* $Id: m-conc.c,v 1.86 2007/05/21 20:16:11 enlf Exp $ */ /* Copyright (C) 1997-2007 by Erik Nolf This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 59 Place - Suite 330, Boston, MA 02111-1307, USA. */ /**************************************************************************** H E A D E R S ****************************************************************************/ #include "m-depend.h" #include #include #include #include #define __USE_XOPEN #include #ifdef HAVE_STDLIB_H #include #endif #ifdef HAVE_STRPTIME #define __USE_XOPEN_EXTENDED #endif #ifdef HAVE_STRING_H #include #endif #ifdef HAVE_STRINGS_H #ifndef _WIN32 #include #endif #endif #ifdef HAVE_UNISTD_H #include #endif #include "medcon.h" /**************************************************************************** D E F I N E S ****************************************************************************/ /* isotope_branching_fraction needed for activity concentration ?*/ #define MDC_ENABLE_ISOTOPE_BRANCHING_FRACTION 1 char * MdcConcFileTypeNames[MDC_CONC_NUM_FILE_TYPES] = { "Unknown", "List mode", "Sinogram", "Normalization", "Attenuation correction", "Image", "Blank", "Unknown/Reserved", "Mu map", }; char * MdcConcAcqModeNames[MDC_CONC_NUM_ACQ_MODES] = { "Unknown", "Blank", "Emission", "Dynamic", "Gated", "Continous bed motion", "Singles transmission", "Windowed Coincidence transmission", "Non-windowed Coincidence transmission" }; char * MdcConcBedMotionNames[MDC_CONC_NUM_BED_MOTIONS] = { "Static or unknown bed motion", "Continous bed motion" }; char * MdcConcTXSrcNames[MDC_CONC_NUM_TX_SRC_TYPES] = { "Unknown TX source type", "TX point source", "TX line source" }; char * MdcConcDataTypeNames[MDC_CONC_NUM_DATA_TYPES] = { "Unknown", "Signed 8 bit", "Signed 16 bit integer, little endian (Intel, DEC)", "Signed 32 bit integer, little endian (Intel, DEC)", "IEEE Float (32 bit), little endian (Intel, DEC)", "IEEE Float (32 bit), big endian (IBM, Sun)", "Signed 16 bit integer, big endian (IBM, Sun)", "Signed 32 bit integer, big endian (IBM, Sun)" }; char * MdcConcOrderModeNames[MDC_CONC_NUM_ORDER_MODES] = { "Element/Axis/View/Ring_diff - view mode", "Element/View/Axis/Ring_Diff - sinogram mode" }; char * MdcConcRebinTypeNames[MDC_CONC_NUM_REBIN_TYPES] = { "Unknown, or no, algorith type", "Full 3D binning (span and ring difference)", "Single-Slice Rebinning", "Fourier Rebinning" }; char * MdcConcReconTypeNames[MDC_CONC_NUM_RECON_TYPES] = { "Unknown, or no, algorithm type", "Filtered Backprojection", "OSEM 2D", "unused", "unused", "unused", "OSEM 3D followed by MAP" }; char * MdcConcOSEM2DTypeNames[MDC_CONC_NUM_OSEM2D_TYPES] = { "Unweighted OSEM2D reconstruction", "Attenuation weighted OSEM2D reconstruction" }; /* deadtime correction applied to the data set */ char * MdcConcDeadCorrTypeNames[MDC_CONC_NUM_DEAD_CORR_TYPES] = { "No deadtime correction applied", "Global estimate based on singles", "CMS estimate based on singles" }; char * MdcConcAttnCorrNames[MDC_CONC_NUM_ATTN_CORR_TYPES] = { "No attenuation applied", "Point source in TX coincidence", "Point source singles based TX", "Segmented point source in TX coincidence", "Segmented point source singles based TX", "Calculated by geometry", "Non-positron source singles based TX" }; char * MdcConcScatterCorrNames[MDC_CONC_NUM_SCATTER_CORR_TYPES] = { "No scatter correction applied", "Fit of emission tail", "Monte Carlo of emission and transmission data", "Direct calculation from analytical formulas" }; char * MdcConcEventTypeNames[MDC_CONC_NUM_EVENT_TYPES] = { "Unknown event type", "Singles", "Prompt events (coincidences)", "Delay events", "Trues" }; char * MdcConcFilterTypeNames[MDC_CONC_NUM_FILTER_TYPES] = { "No filter", "Ramp filter (backprojection) or no filter", "First-order Butterworth window", "Hanning window", "Hamming window", "Parzen window", "Shepp filter", "Second-order Butterworth window", }; char * MdcConcNormTypeNames[MDC_CONC_NUM_NORM_TYPES] = { "No normalization applied", "Point source inversion", "Point source component based", "Cylinder source inversion", "Cylinder source component based" }; char * MdcConcCalibUnitNames[MDC_CONC_NUM_CALIB_UNITS] = { "Unknown calibration units", "nanoCuries/cc", "bequerels/cc" }; char * MdcConcDoseUnitNames[MDC_CONC_NUM_DOSE_UNITS] = { "Unknown dose units", "mCi", "MBq" }; char * MdcConcSubjectOrientationNames[MDC_CONC_NUM_SUBJECT_ORIENTATIONS] = { "Unknown subject orientation", "Feet first, prone", "Head first, prone", "Feet first, supine", "Head first, supine", "Feet first, right", "Head first, right", "Feet first, left", "Head first, left", }; char * MdcConcLengthUnitNames[MDC_CONC_NUM_LENGTH_UNITS] = { "Unknown length units", "millimeters", "centimeters", "inches", }; char * MdcConcWeightUnitNames[MDC_CONC_NUM_WEIGHT_UNITS] = { "Unknown weight units", "grams", "ounces", "kilograms", "pounds", }; char * MdcConcHdrValueNames[MDC_CONC_NUM_HDR_VALUES] = { "version", "model", "institution", "study", "file_name", "file_type", "acquisition_mode", "bed_motion", "total_frames", "isotope", "isotope_half_life", "isotope_branching_fraction", "transaxial_crystals_per_block", "axial_crystals_per_block", "intrinsic_crystal_offset", "transaxial_blocks", "axial_blocks", "transaxial_crystal_pitch", "axial_crystal_pitch", "radius", "radial_fov", "pt_src_radius", "src_radius", "src_cm_per_rev", "tx_src_type", "pt_src_steps_per_rev", "src_steps_per_rev", "default_projections", "default_transaxial_angles", "crystal_thickness", "depth_of_interaction", "transaxial_bin_size", "axial_plane_size", "lld", "uld", "timing_window", "data_type", "data_order", "span", "ring_difference", "number_of_dimensions", "x_dimension", "y_dimension", "z_dimension", "w_dimension", "delta_elements", "x_filter", "y_filter", "z_filter", "histogram_version", "rebinning_type", "rebinning_version", "recon_algorithm", "recon_version", "map_subsets", "map_osem3d_iterations", "map_iterations", "map_beta", "map_blur_type", "map_prior_type", "map_blur_file", "map_pmatrix_file", "osem2d_method", "osem2d_subsets", "osem2d_iterations", "osem2d_em_iterations", "osem2d_map", "osem2d_x_offset", "osem2d_y_offset", "osem2d_zoom", "deadtime_correction_applied", "decay_correction_applied", "normalization_applied", "normalization_filename", "attenuation_applied", "scatter_correction", "scatter_version", "arc_correction_applied", "rotation", "x_offset", "y_offset", "z_offset", "zoom", "pixel_size", "calibration_units", "calibration_factor", "calibration_branching_fraction", "number_of_singles_rates", "investigator", "operator", "study_identifier", "scan_time", "injected_compound", "dose_units", "dose", "injection_time", "injection_decay_correction", "gate_inputs", "gate_bins", "gate_description", "subject_identifier", "subject_genus", "subject_orientation", "subject_length_units", "subject_length", "subject_weight_units", "subject_weight", "subject_phenotype", "study_model", "anesthesia", "analgesia", "other_drugs", "food_access", "water_access", "end_of_header" }; char * MdcConcBlockValueNames[MDC_CONC_NUM_BLOCK_VALUES] = { "frame", "event_type", "gate", "bed", "bed_offset", "ending_bed_offset", "vertical_bed_offset", "data_file_pointer", "frame_start", "frame_duration", "scale_factor", "minimum", "maximum", "deadtime_correction", "decay_correction", "prompts", "delays", "trues", "prompts_rate", "delays_rate", "singles", "end_of_header" }; #define MDC_INPUT_STRING_SIZE 512 #define MDC_CONC_SUPPORTED_VERSION 001.530 #define MDC_MAX_NUM_GARBAGE_LINES 4 /**************************************************************************** internal functions ****************************************************************************/ static MdcConcHdrValue conc_find_next_hdr_line(FILE * hdr_fp, char ** return_line) { char line[MDC_INPUT_STRING_SIZE]; char token[MDC_INPUT_STRING_SIZE]; int conversion_return_value; char done; char valid = MDC_FALSE; MdcConcHdrValue hdr_value = MDC_CONC_HDR_UNKNOWN; MdcConcHdrValue i_value; done = MDC_FALSE; while (!done) { if ( fgets(line, MDC_INPUT_STRING_SIZE, hdr_fp) == NULL) { /* EOF */ done = MDC_TRUE; valid = MDC_FALSE; hdr_value = MDC_CONC_HDR_EOF; *return_line = NULL; } else if (line[0] != '#') { /* skip comment lines */ done = MDC_TRUE; valid = MDC_TRUE; } } if (valid) { conversion_return_value = sscanf(line, "%s ", token); if (conversion_return_value == EOF) hdr_value = MDC_CONC_HDR_EOF; else if (conversion_return_value <= 0) hdr_value = MDC_CONC_HDR_EOF; else { hdr_value = MDC_CONC_HDR_UNKNOWN; *return_line = NULL; for (i_value = 0; i_value < MDC_CONC_NUM_HDR_VALUES; i_value++) { if (strcasecmp(token, MdcConcHdrValueNames[i_value]) == 0) { hdr_value = i_value; i_value = MDC_CONC_NUM_HDR_VALUES; *return_line = (char *)strdup(line); } } if (hdr_value == MDC_CONC_HDR_UNKNOWN) { /* didn't find anything, return the whole line for error msg */ *return_line = (char *)strdup(line); } } } return hdr_value; } static MdcConcBlockValue conc_find_next_block_line(FILE * hdr_fp, char ** return_line) { char line[MDC_INPUT_STRING_SIZE]; char token[MDC_INPUT_STRING_SIZE]; int conversion_return_value; char done; char valid = MDC_FALSE; MdcConcBlockValue block_value = MDC_CONC_BLOCK_UNKNOWN; MdcConcBlockValue i_value; done = MDC_FALSE; while (!done) { if ( fgets(line, MDC_INPUT_STRING_SIZE, hdr_fp) == NULL) { /* read the next line */ done = MDC_TRUE; valid = MDC_FALSE; block_value = MDC_CONC_BLOCK_EOF; *return_line = NULL; } else { if (line[0] != '#') { /* not done if this is a comment line */ done = MDC_TRUE; valid = MDC_TRUE; } } } if (valid) { conversion_return_value = sscanf(line, "%s ", token); if (conversion_return_value == EOF) block_value = MDC_CONC_BLOCK_EOF; else if (conversion_return_value <= 0) block_value = MDC_CONC_BLOCK_EOF; else { block_value = MDC_CONC_BLOCK_UNKNOWN; *return_line = NULL; for (i_value = 0; i_value < MDC_CONC_NUM_BLOCK_VALUES; i_value++) { if (strcasecmp(token, MdcConcBlockValueNames[i_value]) == 0) { block_value = i_value; i_value = MDC_CONC_NUM_BLOCK_VALUES; *return_line = (char *)strdup(line); } } if (block_value == MDC_CONC_BLOCK_UNKNOWN) { /* didn't find anything, return the whole line for error msg */ *return_line = (char *)strdup(line); } } } return block_value; } static float conc_get_float(char * line, int * return_code) { float return_float; *return_code = sscanf(line, "%*s %f", &return_float); if ((*return_code == EOF) || (*return_code <= 0)) return_float = -1.0; return return_float; } static int conc_get_int(char * line, int * return_code) { int return_int; *return_code = sscanf(line, "%*s %d", &return_int); if ((*return_code == EOF) || (*return_code <= 0)) return_int = -1; return return_int; } static void conc_get_int_float(char * line, int * return_code, int * intp, float * floatp) { *return_code = sscanf(line, "%*s %d %f", intp, floatp); if ((*return_code == EOF) || (*return_code <= 0)) { *intp = -1; *floatp= -1.0; } return; } static void conc_get_float_int(char * line, int * return_code, float * floatp, int * intp) { *return_code = sscanf(line, "%*s %f %d", floatp, intp); if ((*return_code == EOF) || (*return_code <= 0)) { *intp = -1; *floatp= -1.0; } return; } static void conc_get_Int32_Int32(char * line, int * return_code, Int32 * num1, Int32 * num2) { *return_code = sscanf(line, "%*s %d %d", num1, num2); if ((*return_code == EOF) || (*return_code <= 0)) *num1 = *num2 = 0; return; } static void conc_get_int_int(char * line, int * return_code, int * num1, int * num2) { *return_code = sscanf(line, "%*s %d %d", num1, num2); if ((*return_code == EOF) || (*return_code <= 0)) *num1 = *num2 = 0; return; } static void conc_get_int_int_float_float(char * line, int * return_code, int * num1, int * num2, float * num3, float * num4) { *return_code = sscanf(line, "%*s %d %d %f %f", num1, num2, num3, num4); if ((*return_code == EOF) || (*return_code <= 0)) *num1 = *num2 = 0; return; } static char * conc_get_string(char * line, int * return_code) { char * return_string; int start_copy; size_t copy_length; /* AML - MacOSX sscan() broken: can't have space before %n, */ /* also added while loop below */ /* *return_code = sscanf(line, "%*s %n", &start_copy) */ *return_code = sscanf(line,"%*s%n",&start_copy); if (*return_code == EOF) return_string = NULL; else { while (line[start_copy] == ' ') start_copy++;/* while added for MacOSX */ copy_length = strcspn(&(line[start_copy]), "\n"); /* return_string = (char *) strndup(&(line[start_copy]), copy_length); */ /* eNlf: strndup replaced for portability */ MdcRemoveEnter(&line[start_copy]); /* remove '\n' as well as '\r */ return_string = malloc(copy_length + 1); if (return_string != NULL) { strncpy(return_string,&(line[start_copy]),copy_length); return_string[copy_length]='\0'; } } return return_string; } static int conc_get_int_string(char * line, int * return_code, int * num1, char ** string1) { int return_int; int conversion_end; *return_code = sscanf(line, "%*s %d%n", &return_int, &conversion_end); if ((*return_code == EOF) || (*return_code <= 0)) return_int = -1; /* the minus 1 is a quick hack to get conc_get_string to work for this */ *string1 = conc_get_string(&(line[conversion_end-1]), return_code); return return_int; } static float conc_convert_injected_dose_to_MBq(float input_id, MdcConcDoseUnits dose_units) { float id; switch(dose_units) { case MDC_CONC_DOSE_UNITS_MILLICURIES: id = MdcmCi2MBq(input_id); break; case MDC_CONC_DOSE_UNITS_UNKNOWN: case MDC_CONC_DOSE_UNITS_MEGA_BEQUERELS: default: id = input_id; break; } return id; } static float conc_convert_weight_to_kg(float input_weight, MdcConcWeightUnits weight_units) { float weight; switch(weight_units) { case MDC_CONC_WEIGHT_UNITS_GRAMS: weight = input_weight/1000.; break; case MDC_CONC_WEIGHT_UNITS_OUNCES: input_weight /= 16.; case MDC_CONC_WEIGHT_UNITS_POUNDS: weight = input_weight / 2.2046226; break; case MDC_CONC_WEIGHT_UNITS_KILOGRAMS: case MDC_CONC_WEIGHT_UNITS_UNKNOWN: default: weight = input_weight; } return weight; } static Int16 conc_save_type(FILEINFO *fi) { Int16 type; /* currently supported BIT8_S, BIT16_S, BIT32_S, FLT32 */ if (MDC_FORCE_INT != MDC_NO) { switch (MDC_FORCE_INT) { case BIT8_U : MdcPrntWarn("CONC Format doesn't support Uint8 type"); case BIT16_S: type = BIT16_S; break; default : type = BIT16_S; } }else{ switch(fi->type) { case BIT8_S : type = BIT8_S; break; case BIT8_U : case BIT16_S: type = BIT16_S; break; case BIT16_U: case BIT32_S: type = BIT32_S; break; case BIT32_U: case BIT64_S: case BIT64_U: case FLT32 : case FLT64 : default : type = FLT32; } } return(type); } /**************************************************************************** F U N C T I O N S ****************************************************************************/ const char *MdcLoadPlaneCONC(FILEINFO *fi, int img) { size_t bytes; IMG_DATA * plane; plane = &fi->image[img]; if (plane->load_location < 0) return("CONC Incorrect plane location in file"); if (plane->buf != NULL) return("CONC Tried to reload plane"); if (fseek(fi->ifp_raw, plane->load_location, SEEK_SET) < 0) { fi->truncated=MDC_YES; return("CONC Could not seek to appropriate file location, truncated read"); } bytes = plane->width*plane->height; bytes *= MdcType2Bytes(plane->type); plane->buf = MdcGetImgBuffer(bytes); if (fread(plane->buf,1,bytes,fi->ifp_raw) != bytes) { fi->truncated=MDC_YES; return("CONC Truncated file read"); } return NULL; } const char *MdcLoadHeaderCONC(FILEINFO *fi) { FILE *hdr_fp = fi->ifp; IMG_DATA * first_plane; IMG_DATA * plane; DYNAMIC_DATA * dd = NULL; BED_DATA * bd = NULL; GATED_DATA * gd = NULL; MdcConcHdrValue hdr_value; MdcConcBlockValue block_value; MdcConcFileTypes file_type=MDC_CONC_FILE_IMAGE; MdcConcWeightUnits weight_units=MDC_CONC_WEIGHT_UNITS_UNKNOWN; MdcConcDoseUnits dose_units = MDC_CONC_DOSE_UNITS_UNKNOWN; MdcConcAcqModes acq_type=MDC_CONC_ACQ_UNKNOWN; char * line = NULL; char done; float temp_float, temp_float2; int temp_int, temp_int2; char * temp_string; char * raw_filename=NULL; char base_filename[MDC_MAX_PATH+1]; char * header_derived_filename=NULL; char * pfilename; int return_code; MdcConcDeadCorrTypes deadtime_correction; MdcConcCalibUnits calibration_units=0; MdcConcReconTypes recon_type=-1; float osem2d_recon_zoom=-1; float calibration_factor=1.0; float isotope_branching_factor=1.0; int i_bed, i_gate, i_frame, i_plane, img; Int32 high_file_pointer, low_file_pointer; Uint32 number; char found_frames = MDC_FALSE; char found_data_type = MDC_FALSE; struct tm scan_time; /* MARK: unused yet #ifdef LONG_8BYTE Uint64 location; #endif */ int num_garbage_lines = 0; char failed_to_read_64_bit=MDC_FALSE; long plane_bytes; char transmission_scan = MDC_FALSE; if (MDC_VERBOSE) MdcPrntMesg("CONC Reading <%s> ...",fi->ifname); /* initialize */ fi->modality = M_PT; /* read through the header, looking at all the tokens */ done = MDC_FALSE; while (!done) { hdr_value = conc_find_next_hdr_line(hdr_fp, &line); switch (hdr_value) { case MDC_CONC_HDR_VERSION: if (MDC_INFO) MdcPrntScrn("Concorde file version:\t\t%f\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_MODEL: if (MDC_INFO) MdcPrntScrn("Scanner model:\t\t\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_INSTITUTION: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Institution:\t\t\t%s\n",temp_string); if (strlen(temp_string)) MdcStringCopy(fi->institution,temp_string,strlen(temp_string)); MdcFree(temp_string); break; case MDC_CONC_HDR_STUDY: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Study:\t\t\t\t%s\n",temp_string); MdcStringCopy(fi->study_id, temp_string,strlen(temp_string)); MdcFree(temp_string); break; case MDC_CONC_HDR_FILE_NAME: raw_filename = conc_get_string(line, &return_code); if (MDC_INFO) MdcPrntScrn("Raw data file:\t\t\t%s\n",raw_filename); pfilename = strrchr(raw_filename,'/'); if (pfilename == NULL) pfilename = strrchr(raw_filename,'\\'); if (pfilename != NULL) pfilename++; else pfilename = raw_filename; base_filename[0]='\0'; if (fi->idir != NULL) { strncpy(base_filename, fi->idir, MDC_MAX_PATH); strncat(base_filename, MDC_PATH_DELIM_STR, MDC_MAX_PATH); } strncat(base_filename, pfilename, MDC_MAX_PATH); if (MDC_INFO) MdcPrntScrn("Base name of raw data file:\t%s\n",pfilename); break; case MDC_CONC_HDR_FILE_TYPE: temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_FILE_TYPES)) { file_type = MDC_CONC_FILE_UNKNOWN; }else{ file_type = temp_int; } if (MDC_INFO) MdcPrntScrn("File type:\t\t\t%d=%s\n", temp_int,MdcConcFileTypeNames[file_type]); switch(file_type) { case MDC_CONC_FILE_MU_MAP: case MDC_CONC_FILE_IMAGE: fi->reconstructed = MDC_YES; break; case MDC_CONC_FILE_ATTENUATION: case MDC_CONC_FILE_SINOGRAM: case MDC_CONC_FILE_NORMALIZATION: fi->reconstructed = MDC_NO; break; default: return("CONC Cannot handle this Concorde/uPET file type"); break; } break; case MDC_CONC_HDR_ACQUISITION_MODE: temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_ACQ_MODES)) { acq_type = MDC_CONC_ACQ_UNKNOWN; }else{ acq_type = temp_int; } if (MDC_INFO) MdcPrntScrn("Acquisition type:\t\t%d=%s\n", temp_int,MdcConcAcqModeNames[acq_type]); switch(acq_type) { case MDC_CONC_ACQ_BLANK: fi->acquisition_type = MDC_ACQUISITION_STATIC; break; case MDC_CONC_ACQ_EMISSION: fi->acquisition_type = MDC_ACQUISITION_TOMO; break; case MDC_CONC_ACQ_DYNAMIC: fi->acquisition_type = MDC_ACQUISITION_DYNAMIC; break; case MDC_CONC_ACQ_GATED: fi->acquisition_type = MDC_ACQUISITION_GATED; break; case MDC_CONC_ACQ_CONTINUOUS: fi->acquisition_type = MDC_ACQUISITION_TOMO; break; case MDC_CONC_ACQ_SINGLES: fi->acquisition_type = MDC_ACQUISITION_STATIC; transmission_scan = MDC_TRUE; break; case MDC_CONC_ACQ_WINDOWED_COINCIDENCE: case MDC_CONC_ACQ_NON_WINDOWED_COINCIDENCE: fi->acquisition_type = MDC_ACQUISITION_TOMO; transmission_scan = MDC_TRUE; break; case MDC_CONC_ACQ_UNKNOWN: default: fi->acquisition_type = MDC_ACQUISITION_UNKNOWN; break; } break; case MDC_CONC_HDR_BED_MOTION: { MdcConcBedMotion bed_motion; temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_BED_MOTIONS)) { bed_motion = MDC_CONC_BED_MOTION_STATIC; }else{ bed_motion = temp_int; } if (MDC_INFO) MdcPrntScrn("Bed Motion type:\t\t%d=%s\n", temp_int,MdcConcBedMotionNames[bed_motion]); if (bed_motion != MDC_CONC_BED_MOTION_STATIC) { MdcPrntWarn("CONC Don't know how to handle bed motion:\t%s", MdcConcBedMotionNames[bed_motion]); } } break; /* don't care */ case MDC_CONC_HDR_TOTAL_FRAMES: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Number of frames:\t\t%d\n",temp_int); if (temp_int <= 0) return("CONC Header reported no frames of data"); fi->dim[4] = temp_int; found_frames = MDC_TRUE; break; case MDC_CONC_HDR_ISOTOPE: temp_string = conc_get_string(line, &return_code); if (MDC_INFO) MdcPrntScrn("Isotope used:\t\t\t%s\n",temp_string); MdcStringCopy(fi->isotope_code, temp_string, strlen(temp_string)); MdcFree(temp_string); break; case MDC_CONC_HDR_ISOTOPE_HALF_LIFE: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Isotope half life:\t\t%5.3f\n",temp_float); fi->isotope_halflife = temp_float; break; case MDC_CONC_HDR_ISOTOPE_BRANCHING_FRACTION: temp_float = conc_get_float(line, &return_code); isotope_branching_factor = temp_float; if (MDC_INFO) MdcPrntScrn("Isotope branching fraction:\t%5.3f\n",temp_float); break; case MDC_CONC_HDR_TRANSAXIAL_CRYSTALS_PER_BLOCK: if (MDC_INFO) MdcPrntScrn("Transaxial crystals per block:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_AXIAL_CRYSTALS_PER_BLOCK: if (MDC_INFO) MdcPrntScrn("Axial crystals per block:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_INTRINSIC_CRYSTAL_OFFSET: if (MDC_INFO) MdcPrntScrn("Crystal off. for intrinsic rot:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_TRANSAXIAL_BLOCKS: if (MDC_INFO) MdcPrntScrn("Transaxial blocks:\t\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_AXIAL_BLOCKS: if (MDC_INFO) MdcPrntScrn("Axial blocks:\t\t\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_TRANSAXIAL_CRYSTAL_PITCH: temp_float = 10.0*conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Transaxial crystal pitch:\t%5.3f (mm)\n",temp_float); break; case MDC_CONC_HDR_AXIAL_CRYSTAL_PITCH: temp_float = 10.0*conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Axial crystal pitch:\t\t%5.3f (mm)\n",temp_float); if (MDC_INFO) MdcPrntScrn("Slice width (z):\t\t%5.3f (mm)\n",temp_float/2.0); fi->pixdim[3] = temp_float/2.0; break; case MDC_CONC_HDR_RADIUS: if (MDC_INFO) MdcPrntScrn("Ring radius (to crystal face):\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_RADIAL_FOV: if (MDC_INFO) MdcPrntScrn("Radial field of view:\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_PT_SRC_RADIUS: case MDC_CONC_HDR_SRC_RADIUS: if (MDC_INFO) MdcPrntScrn("Source radius:\t\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_SRC_CM_PER_REV: if (MDC_INFO) MdcPrntScrn("Source CM/Rev:\t\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_TX_SRC_TYPE: if ((acq_type == MDC_CONC_ACQ_SINGLES) || (acq_type == MDC_CONC_ACQ_WINDOWED_COINCIDENCE) || (acq_type == MDC_CONC_ACQ_NON_WINDOWED_COINCIDENCE)) { MdcConcTXSrcTypes tx_src; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_TX_SRC_TYPES)) { tx_src = temp_int; }else{ tx_src = MDC_CONC_TX_SRC_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Transmission source:\t\t%d=%s\n", temp_int, MdcConcTXSrcNames[tx_src]); } break; case MDC_CONC_HDR_PT_SRC_STEPS_PER_REV: case MDC_CONC_HDR_SRC_STEPS_PER_REV: if (MDC_INFO) MdcPrntScrn("Source encoder steps/rev.:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_DEFAULT_PROJECTIONS: if (MDC_INFO) MdcPrntScrn("Default # of projections:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_DEFAULT_TRANSAXIAL_ANGLES: if (MDC_INFO) MdcPrntScrn("Default # of transaxial angles:\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_CRYSTAL_THICKNESS: if (MDC_INFO) MdcPrntScrn("Crystal thickness:\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_DEPTH_OF_INTERACTION: if (MDC_INFO) MdcPrntScrn("Depth of interaction:\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_TRANSAXIAL_BIN_SIZE: if (MDC_INFO) MdcPrntScrn("Transaxial bin size:\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_AXIAL_PLANE_SIZE: if (MDC_INFO) MdcPrntScrn("Axial plane size:\t\t%5.3f (mm)\n", 10.0*conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_LLD: if (MDC_INFO) MdcPrntScrn("Lower level energy threshold:\t%5.3f (KeV)\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_ULD: if (MDC_INFO) MdcPrntScrn("Upper level energy threshold:\t%5.3f (KeV)\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_TIMING_WINDOW: if (MDC_INFO) MdcPrntScrn("Coincidence timing window:\t%d (nsecs)\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_DATA_TYPE: { MdcConcDataTypes data_type; temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_DATA_TYPES)) { data_type = MDC_CONC_DATA_UNKNOWN; }else{ data_type = temp_int; } if (MDC_INFO) MdcPrntScrn("Data type:\t\t\t%d=%s\n", temp_int,MdcConcDataTypeNames[data_type]); fi->endian=MDC_FILE_ENDIAN= MDC_LITTLE_ENDIAN; switch(data_type) { case MDC_CONC_DATA_SBYTE: fi->bits = 8; fi->type = BIT8_S; break; case MDC_CONC_DATA_SSHORT_BE: fi->endian=MDC_FILE_ENDIAN = MDC_BIG_ENDIAN; case MDC_CONC_DATA_SSHORT_LE: fi->bits = 16; fi->type = BIT16_S; break; case MDC_CONC_DATA_SINT_BE: fi->endian=MDC_FILE_ENDIAN = MDC_BIG_ENDIAN; case MDC_CONC_DATA_SINT_LE: fi->bits = 32; fi->type = BIT32_S; break; case MDC_CONC_DATA_FLOAT_BE: fi->endian=MDC_FILE_ENDIAN = MDC_BIG_ENDIAN; case MDC_CONC_DATA_FLOAT_LE: fi->bits = 32; fi->type = FLT32; break; case MDC_CONC_DATA_UNKNOWN: case MDC_CONC_NUM_DATA_TYPES: default: return("CONC Unknown data type for Concorde/uPET file"); } found_data_type = MDC_TRUE; } break; case MDC_CONC_HDR_DATA_ORDER: { MdcConcOrderModes mode_type; temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_ORDER_MODES)) { mode_type = MDC_CONC_ORDER_SINOGRAM; }else{ mode_type = temp_int; } if (MDC_INFO) MdcPrntScrn("Data order:\t\t\t%d=%s\n", temp_int,MdcConcOrderModeNames[mode_type]); if (mode_type != MDC_CONC_ORDER_SINOGRAM) return("CONC Can only handle Concorde/uPET scans in x/y/z/w order"); } break; case MDC_CONC_HDR_SPAN: if (MDC_INFO) MdcPrntScrn("Span of data set:\t\t%d\n", conc_get_int(line, &return_code)); break; case MDC_CONC_HDR_RING_DIFFERENCE: if (MDC_INFO) MdcPrntScrn("Maximum ring difference:\t%d\n", conc_get_int(line, &return_code)); break; /* don't care */ case MDC_CONC_HDR_NUMBER_OF_DIMENSIONS: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Number of dimensions:\t\t%d\n",temp_int); /* Concorde doesn't use this consistently, so ignore it */ break; case MDC_CONC_HDR_X_DIMENSION: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("X dimension:\t\t\t%d\n",temp_int); fi->dim[1]=temp_int; break; case MDC_CONC_HDR_Y_DIMENSION: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Y dimension:\t\t\t%d\n",temp_int); fi->dim[2]=temp_int; break; case MDC_CONC_HDR_Z_DIMENSION: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Z dimension:\t\t\t%d\n",temp_int); if ((file_type == MDC_CONC_FILE_SINOGRAM) || (file_type == MDC_CONC_FILE_ATTENUATION) || (file_type == MDC_CONC_FILE_NORMALIZATION)) fi->dim[3]=0; else fi->dim[3]=temp_int; break; case MDC_CONC_HDR_W_DIMENSION: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("W dimension (ignored):\t\t%d\n",temp_int); break; case MDC_CONC_HDR_DELTA_ELEMENTS: { int delta, planes; conc_get_int_int(line, &return_code, &delta, &planes); if (MDC_INFO) MdcPrntScrn("\tdelta elements (%d):\t%d\n", delta, planes); fi->dim[3] += planes; } break; case MDC_CONC_HDR_X_FILTER: { MdcConcFilterTypes filter_type; conc_get_int_float(line, &return_code, &temp_int, &temp_float); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_FILTER_TYPES)) { filter_type = temp_int; }else{ filter_type = MDC_CONC_FILTER_NONE; } if (MDC_INFO) MdcPrntScrn("Filter (X):\t\t\t%d=%s\n", temp_int, MdcConcFilterTypeNames[filter_type]); MdcStringCopy(fi->filter_type, MdcConcFilterTypeNames[filter_type], strlen(MdcConcFilterTypeNames[filter_type])); } break; case MDC_CONC_HDR_Y_FILTER: break; /* don't care */ case MDC_CONC_HDR_Z_FILTER: break; /* don't care */ case MDC_CONC_HDR_HISTOGRAM_VERSION: break; /* don't care */ case MDC_CONC_HDR_REBINNING_TYPE: break; /* don't care */ case MDC_CONC_HDR_REBINNING_VERSION: break; /* don't care */ case MDC_CONC_HDR_RECON_ALGORITHM: temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_RECON_TYPES)) { recon_type = temp_int; }else{ recon_type = MDC_CONC_RECON_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Reconstruction type:\t\t%d=%s\n", temp_int, MdcConcReconTypeNames[recon_type]); MdcStringCopy(fi->recon_method, MdcConcReconTypeNames[recon_type], strlen(MdcConcReconTypeNames[recon_type])); break; case MDC_CONC_HDR_RECON_VERSION: if (MDC_INFO) MdcPrntScrn("Reconstruction version:\t\t%f\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_MAP_SUBSETS: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM3D Subsets:\t\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_MAP_OSEM3D_ITERATIONS: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM3D Iterations:\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_MAP_ITERATIONS: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP Iterations:\t\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_MAP_BETA: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP Beta:\t\t\t%5.3f\n", temp_float); } break; case MDC_CONC_HDR_MAP_BLUR_TYPE: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP Blur Type:\t\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_MAP_PRIOR_TYPE: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP Prior Type:\t\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_MAP_BLUR_FILE: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_string = conc_get_string(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP Blur File:\t\t\t%s\n",temp_string); MdcFree(temp_string); } break; case MDC_CONC_HDR_MAP_PMATRIX_FILE: if (recon_type == MDC_CONC_RECON_OSEM3D_MAP) { temp_string = conc_get_string(line, &return_code); if (MDC_INFO) MdcPrntScrn("MAP PMatrix File:\t\t%s\n",temp_string); MdcFree(temp_string); } break; case MDC_CONC_HDR_OSEM2D_METHOD: if (recon_type == MDC_CONC_RECON_OSEM2D) { MdcConcOSEM2DTypes osem_type; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_OSEM2D_TYPES)) { osem_type = temp_int; if (MDC_INFO) MdcPrntScrn("OSEM2D type:\t\t\t%d=%s\n", temp_int, MdcConcOSEM2DTypeNames[osem_type]); }else{ osem_type = MDC_CONC_OSEM2D_UNKNOWN; } } break; case MDC_CONC_HDR_OSEM2D_SUBSETS: if (recon_type == MDC_CONC_RECON_OSEM2D) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM2D Subsets:\t\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_OSEM2D_ITERATIONS: if (recon_type == MDC_CONC_RECON_OSEM2D) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM2D Iterations:\t\t%d\n", temp_int); } break; case MDC_CONC_HDR_OSEM2D_EM_ITERATIONS: if (recon_type == MDC_CONC_RECON_OSEM2D) { temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM2D EM Iterations:\t\t%d\n", temp_int); } break; /* Epsilon and power values for map regularization (float integer) */ case MDC_CONC_HDR_OSEM2D_MAP: if (recon_type == MDC_CONC_RECON_OSEM2D) { conc_get_float_int(line, &return_code, &temp_float, &temp_int); temp_int = conc_get_int(line, &return_code); if (MDC_INFO) { MdcPrntScrn("OSEM2D Epsilon:\t\t\t%5.3f\n", temp_float); MdcPrntScrn("OSEM2D Power:\t\t\t%d\n", temp_int); } } break; case MDC_CONC_HDR_OSEM2D_X_OFFSET: if (recon_type == MDC_CONC_RECON_OSEM2D) { temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM2D X_Offset:\t\t%5.3f\n", temp_float); } break; case MDC_CONC_HDR_OSEM2D_Y_OFFSET: if (recon_type == MDC_CONC_RECON_OSEM2D) { temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("OSEM2D Y_Offset:\t\t%5.3f\n", temp_float); } break; case MDC_CONC_HDR_OSEM2D_ZOOM: osem2d_recon_zoom = conc_get_float(line, &return_code); if (recon_type == MDC_CONC_RECON_OSEM2D) { if (MDC_INFO) MdcPrntScrn("OSEM2D Recon Zoom:\t\t%5.3f\n", osem2d_recon_zoom); } break; case MDC_CONC_HDR_DEADTIME_CORRECTION_APPLIED: temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_DEAD_CORR_TYPES)) { deadtime_correction = temp_int; }else{ deadtime_correction = MDC_CONC_DEAD_CORR_NONE; } if (MDC_INFO) MdcPrntScrn("Deadtime correction applied:\t%d=%s\n", temp_int, MdcConcDeadCorrTypeNames[deadtime_correction]); break; case MDC_CONC_HDR_DECAY_CORRECTION_APPLIED: fi->decay_corrected = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Decay correction applied:\t%s\n", fi->decay_corrected ? "applied" : "not applied"); break; case MDC_CONC_HDR_NORMALIZATION_APPLIED: { MdcConcNormTypes norm_type; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_NORM_TYPES)) { norm_type = temp_int; }else{ norm_type = MDC_CONC_NORM_NONE; } if (MDC_INFO) MdcPrntScrn("Normalization applied:\t\t%d=%s\n", temp_int, MdcConcNormTypeNames[norm_type]); } break; case MDC_CONC_HDR_NORMALIZATION_FILENAME: temp_string = conc_get_string(line, &return_code); if (MDC_INFO) MdcPrntScrn("Normalization file:\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_ATTENUATION_APPLIED: { MdcConcAttnCorrTypes attn_type; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_ATTN_CORR_TYPES)) { attn_type = temp_int; }else{ attn_type = MDC_CONC_ATTN_CORR_NONE; } if (MDC_INFO) MdcPrntScrn("Attenuation correction applied:\t%d=%s\n", temp_int, MdcConcAttnCorrNames[attn_type]); } break; case MDC_CONC_HDR_SCATTER_CORRECTION: { MdcConcScatterCorrTypes scatter_type; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_SCATTER_CORR_TYPES)) { scatter_type = temp_int; }else{ scatter_type = MDC_CONC_SCATTER_CORR_NONE; } if (MDC_INFO) MdcPrntScrn("Scatter correction applied:\t%d=%s\n", temp_int, MdcConcScatterCorrNames[scatter_type]); } break; case MDC_CONC_HDR_SCATTER_VERSION: if (MDC_INFO) MdcPrntScrn("Scatter Correction version:\t%f\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_ARC_CORRECTION_APPLIED: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Arc correction applied:\t\t%s\n", temp_int ? "applied" : "not applied"); break; case MDC_CONC_HDR_ROTATION: break; /* don't care */ case MDC_CONC_HDR_X_OFFSET: break; /* don't care */ case MDC_CONC_HDR_Y_OFFSET: break; /* don't care */ case MDC_CONC_HDR_Z_OFFSET: break; /* don't care */ case MDC_CONC_HDR_ZOOM: break; /* don't care */ case MDC_CONC_HDR_PIXEL_SIZE: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Pixel Size (x,y) (mm):\t\t%5.3f\n", 10*temp_float); fi->pixdim[1] = fi->pixdim[2] = 10*temp_float; break; case MDC_CONC_HDR_CALIBRATION_UNITS: temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_CALIB_UNITS)) { calibration_units = temp_int; }else{ calibration_units = MDC_CONC_CALIB_UNITS_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Calibration Units:\t\t%d=%s\n", temp_int, MdcConcCalibUnitNames[calibration_units]); break; case MDC_CONC_HDR_CALIBRATION_FACTOR: temp_float = conc_get_float(line, &return_code); calibration_factor *= temp_float; if (MDC_INFO) MdcPrntScrn("Calibration Factor:\t\t%5.3f\n",temp_float); break; case MDC_CONC_HDR_CALIBRATION_BRANCHING_FRACTION: if (MDC_INFO) MdcPrntScrn("Calibration branching fraction:\t%5.3f\n", conc_get_float(line, &return_code)); break; case MDC_CONC_HDR_NUMBER_OF_SINGLES_RATES: break; /* don't care */ case MDC_CONC_HDR_INVESTIGATOR: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Investigator:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_OPERATOR: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Operator:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_STUDY_IDENTIFIER: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Study ID:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_SCAN_TIME: #ifdef HAVE_STRPTIME temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Scan Time:\t\t\t%s\n",temp_string); if (strlen(temp_string)) { scan_time.tm_mday=0; scan_time.tm_mon=0; scan_time.tm_year=0; scan_time.tm_hour=0; scan_time.tm_min=0; scan_time.tm_sec=0; strptime(temp_string, "%a %b %d %T %Y", &scan_time); if (scan_time.tm_year != 0) { fi->study_date_day = scan_time.tm_mday; fi->study_date_month = scan_time.tm_mon+1; fi->study_date_year = scan_time.tm_year+1900; fi->study_time_hour = scan_time.tm_hour; fi->study_time_minute = scan_time.tm_min; fi->study_time_second = scan_time.tm_sec; } } MdcFree(temp_string); #endif break; case MDC_CONC_HDR_INJECTED_COMPOUND: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Injected Compound:\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_DOSE_UNITS: dose_units = conc_get_int(line, &return_code); if ((dose_units < 0) || (dose_units >= MDC_CONC_NUM_DOSE_UNITS)) { dose_units = MDC_CONC_DOSE_UNITS_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Dose Units:\t\t\t%d=%s\n", dose_units, MdcConcDoseUnitNames[dose_units]); /* next line incase INJECTED_DOSE read before DOSE_UNITS */ fi->injected_dose= conc_convert_injected_dose_to_MBq(fi->injected_dose, dose_units); break; case MDC_CONC_HDR_INJECTED_DOSE: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Injected Dose:\t\t\t%5.3f\n", temp_float); fi->injected_dose= conc_convert_injected_dose_to_MBq(temp_float,dose_units); break; case MDC_CONC_HDR_INJECTION_TIME: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Injection Time:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_INJECTION_DECAY_CORRECTION: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Injection Decay Correction:\t%5.3f\n",temp_float); break; case MDC_CONC_HDR_GATE_INPUTS: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("Number of Gate Inputs:\t\t%d\n", temp_int); MdcGetStructGD(fi, (unsigned)temp_int); break; case MDC_CONC_HDR_GATE_BINS: conc_get_int_int_float_float(line, &return_code, &temp_int, &temp_int2, &temp_float, &temp_float2); if (MDC_INFO) { MdcPrntScrn("Number of Bins for Gate %d:\t%d\n", temp_int, temp_int2); MdcPrntScrn(" min/max gate cycle:\t%5.3f/%5.3f\n", temp_float, temp_float2); } if (fi->gdata == NULL) return("CONC gate_inputs line must preceed gate_bins line"); if (temp_int >= fi->gatednr) return("CONC more gates found then specified in gate_inputs line"); gd = &(fi->gdata[temp_int]); gd->nr_projections = temp_int2; gd->window_low = temp_float*1000.; gd->window_high = temp_float*1000.; if (fi->dim[5] == 0) fi->dim[5] = temp_int2; else fi->dim[5]*=temp_int2; fi->dim[4]/=temp_int2; /* Concorde sets frames = total gates */ break; case MDC_CONC_HDR_GATE_DESCRIPTION: conc_get_int_string(line, &return_code, &temp_int, &temp_string); if (MDC_INFO) MdcPrntScrn("Gate %d description:\t\t%s\n", temp_int, temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_SUBJECT_IDENTIFIER: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Subject ID:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_SUBJECT_GENUS: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Subject Genus:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_SUBJECT_ORIENTATION: temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_SUBJECT_ORIENTATIONS)) { temp_int = MDC_CONC_SUBJECT_ORIENTATION_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Subject Orientation:\t\t%d=%s\n", temp_int, MdcConcSubjectOrientationNames[temp_int]); switch(temp_int) { case MDC_CONC_SUBJECT_ORIENTATION_FEET_PRONE: fi->pat_slice_orient=MDC_PRONE_FEETFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_HEAD_PRONE: fi->pat_slice_orient=MDC_PRONE_HEADFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_FEET_SUPINE: fi->pat_slice_orient=MDC_SUPINE_FEETFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_HEAD_SUPINE: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_SAGITTAL; break; case MDC_CONC_SUBJECT_ORIENTATION_FEET_RIGHT: fi->pat_slice_orient=MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_HEAD_RIGHT: fi->pat_slice_orient=MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_FEET_LEFT: fi->pat_slice_orient=MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_HEAD_LEFT: fi->pat_slice_orient=MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL; break; case MDC_CONC_SUBJECT_ORIENTATION_UNKNOWN: default: fi->pat_slice_orient=MDC_UNKNOWN; break; } break; case MDC_CONC_HDR_SUBJECT_LENGTH_UNITS: temp_int = conc_get_int(line, &return_code); if ((temp_int < 0) || (temp_int >= MDC_CONC_NUM_LENGTH_UNITS)) { temp_int = MDC_CONC_LENGTH_UNITS_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Length Units:\t\t\t%d=%s\n", temp_int, MdcConcLengthUnitNames[temp_int]); break; case MDC_CONC_HDR_SUBJECT_LENGTH: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Subject Length:\t\t\t%5.3f\n", temp_float); break; case MDC_CONC_HDR_SUBJECT_WEIGHT_UNITS: weight_units = conc_get_int(line, &return_code); if ((weight_units < 0) || (weight_units >= MDC_CONC_NUM_WEIGHT_UNITS)) { weight_units = MDC_CONC_WEIGHT_UNITS_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("Weight Units:\t\t\t%d=%s\n", weight_units, MdcConcWeightUnitNames[weight_units]); /* next line incase SUBJECT_WEIGHT read before WEIGHT_UNITS */ fi->patient_weight = conc_convert_weight_to_kg(fi->patient_weight, weight_units); break; case MDC_CONC_HDR_SUBJECT_WEIGHT: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("Subject Weight:\t\t\t%5.3f\n", temp_float); fi->patient_weight = conc_convert_weight_to_kg(temp_float,weight_units); break; case MDC_CONC_HDR_SUBJECT_PHENOTYPE: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Subject Phenotype:\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_STUDY_MODEL: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Study Model:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_ANESTHESIA: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Anesthesia:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_ANALGESIA: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Analgesia:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_OTHER_DRUGS: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Other drugs:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_FOOD_ACCESS: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Food Access:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_WATER_ACCESS: temp_string = conc_get_string(line, &return_code); if ((MDC_INFO) && (strlen(temp_string))) MdcPrntScrn("Water Access:\t\t\t%s\n",temp_string); MdcFree(temp_string); break; case MDC_CONC_HDR_END_OF_HEADER: done = MDC_TRUE; break; case MDC_CONC_HDR_EOF: done = MDC_TRUE; return("CONC Got inappropriate EOF on reading Concorde/uPET header"); break; case MDC_CONC_HDR_UNKNOWN: default: if (num_garbage_lines < MDC_MAX_NUM_GARBAGE_LINES) { MdcPrntWarn("CONC Uninterpretable line: %s",line); } num_garbage_lines++; break; } MdcFree(line); } if (!found_frames || !found_data_type) { return("CONC Header file didn't contain necessary entries"); } if (MDC_ECHO_ALIAS == MDC_YES) { MdcEchoAliasName(fi); return(NULL); } /* figure out remaining parameters of the FILEINFO structure */ fi->dim[0]=5; fi->pixdim[0]=5; fi->dim[6]=1; /* we can't do multiple bed positions at the moment */ fi->dim[7]=1; number = fi->dim[7]*fi->dim[6]*fi->dim[5]*fi->dim[4]*fi->dim[3]; fi->truncated = MDC_NO; if (number == 0) return("CONC No valid images specified"); /* not worrying about filling in patient orientation */ /* not worrying about patient sex */ /* figure out a guess at the raw filename from the header name */ /* allows us to complain if the raw filename in the header is different */ MdcMergePath(fi->ipath,fi->idir,fi->ifname); header_derived_filename = (char *)strdup(fi->ipath); MdcSplitPath(fi->ipath,fi->idir,fi->ifname); pfilename = strstr(header_derived_filename, ".hdr"); if (pfilename == NULL) { MdcFree(header_derived_filename); header_derived_filename = NULL; } else { strcpy(pfilename, "\0"); } /* complain if the filename derived from the header's filename doesn't match the file_name entry in the header */ if (header_derived_filename != NULL) { if (strcmp(header_derived_filename, base_filename) != 0) { MdcPrntWarn("CONC Name mismatch between header entry & derived file\n" \ "\t\t using: %s\n", base_filename); } } #if MDC_ENABLE_ISOTOPE_BRANCHING_FRACTION /* adjust the calibration factor with the branching factor if needed */ if (!transmission_scan) calibration_factor /= isotope_branching_factor; #endif /* need to open up the raw data file */ fi->ifp_raw = fopen(raw_filename, "rb"); /* try filename with bogus path info removed, *nix or DOS style */ if (fi->ifp_raw == NULL) fi->ifp_raw = fopen(base_filename, "rb"); /* maybe a compressed file */ if (fi->ifp_raw == NULL) { MdcAddCompressionExt(fi->compression, base_filename); if (MdcFileExists(base_filename)) { if (MdcDecompressFile(base_filename) != MDC_OK) { MdcFree(raw_filename); return("CONC Decompression image file failed"); } fi->ifp_raw = fopen(base_filename, "rb"); if (fi->ifp_raw != NULL) unlink(base_filename); /* delete after use */ } } /* allright, try guessing the filename */ if ((fi->ifp_raw == NULL) && (header_derived_filename != NULL)) { fi->ifp_raw = fopen(header_derived_filename, "rb"); if (fi->ifp_raw != NULL) { /* complain we're picking our own data file */ MdcPrntWarn("CONC Header specified raw file (%s) not found," \ "\t\t using: %s",base_filename, header_derived_filename); } } /* fix original path and free strings */ MdcFree(raw_filename); if (header_derived_filename != NULL) MdcFree(header_derived_filename); if (fi->ifp_raw == NULL) return("CONC Couldn't open raw data file"); /* malloc IMG_DATA structs */ if (!MdcGetStructID(fi,number)) return("CONC Bad malloc IMG_DATA structs"); /* malloc DYNAMIC_DATA structs */ if (!MdcGetStructDD(fi,(unsigned) (fi->dim[4]*fi->dim[5]))) return("CONC Bad malloc DYNAMIC_DATA structs"); /* malloc BED_DATA structs */ if (!MdcGetStructBD(fi,(unsigned)fi->dim[6])) return("CONC Bad malloc BED_DATA structs"); img = 0; /* don't yet have support for multiple bed files */ i_bed = 0; /* read each of the gates and frames */ for (i_gate = 0; i_gate < fi->dim[5]; i_gate++) { for (i_frame = 0; i_frame < fi->dim[4]; i_frame++) { if (fi->dynnr > 0) dd = &fi->dyndata[i_frame+i_gate*fi->dim[4]]; if (dd != NULL) dd->nr_of_slices = fi->dim[3]; if (fi->bednr > 0) bd = &fi->beddata[i_bed]; if (MDC_INFO) MdcPrintLine('-', MDC_HALF_LENGTH); if (MDC_INFO) MdcPrntScrn("Gate: \t\t\t\t%d\n",i_gate); if (MDC_INFO) MdcPrntScrn("Frame:\t\t\t\t%d\n",i_frame); first_plane = &fi->image[img]; first_plane->width = fi->dim[1]; first_plane->height = fi->dim[2]; first_plane->bits = fi->bits; first_plane->type = fi->type; first_plane->pixel_xsize = fi->pixdim[1]; first_plane->pixel_ysize = fi->pixdim[2]; first_plane->slice_width = fi->pixdim[3]; if (recon_type == MDC_CONC_RECON_OSEM2D) if (osem2d_recon_zoom > 0) first_plane->recon_scale = osem2d_recon_zoom; /* otherwise use default */ /* continue reading through the header to get the frame info */ done = MDC_FALSE; high_file_pointer=low_file_pointer=-1; while (!done) { block_value = conc_find_next_block_line(hdr_fp, &line); switch (block_value) { case MDC_CONC_BLOCK_FRAME: temp_int = conc_get_int(line, &return_code); if (temp_int != i_frame) { MdcPrntWarn("CONC Detected frame numbering discrepancy"); } break; case MDC_CONC_BLOCK_EVENT_TYPE: { MdcConcEventTypes event_type; temp_int = conc_get_int(line, &return_code); if ((temp_int >= 0) && (temp_int < MDC_CONC_NUM_EVENT_TYPES)) { event_type = temp_int; }else{ event_type = MDC_CONC_EVENT_UNKNOWN; } if (MDC_INFO) MdcPrntScrn("\tEvent type for block:\t%d=%s\n", temp_int, MdcConcEventTypeNames[event_type]); } break; case MDC_CONC_BLOCK_GATE: temp_int = conc_get_int(line, &return_code); if (temp_int != i_gate) { MdcPrntWarn("CONC Detected gate numbering discrepancy"); } break; case MDC_CONC_BLOCK_BED: temp_int = conc_get_int(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tBed #:\t\t\t%d\n",temp_int); if (temp_int != 0) return("CONC Can't handle multiple bed files"); break; case MDC_CONC_BLOCK_BED_OFFSET: temp_float = 10*conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tBed offset (mm):\t%5.3f\n", temp_float); bd->hoffset = temp_float; break; case MDC_CONC_BLOCK_ENDING_BED_OFFSET: if (MDC_INFO) MdcPrntScrn("\tBed ending offset (mm):\t%5.3f\n" ,10*conc_get_float(line, &return_code)); break; case MDC_CONC_BLOCK_VERTICAL_BED_OFFSET: temp_float = 10*conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tBed vert. offset (mm):\t%5.3f\n", temp_float); bd->voffset = temp_float; break; case MDC_CONC_BLOCK_DATA_FILE_POINTER: conc_get_Int32_Int32(line, &return_code, &high_file_pointer , &low_file_pointer); break; case MDC_CONC_BLOCK_FRAME_START: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tFrame start (s):\t%5.3f\n",temp_float); if (dd != NULL) dd->time_frame_start = temp_float * 1000.; break; case MDC_CONC_BLOCK_FRAME_DURATION: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tFrame duration (s):\t%5.3f\n",temp_float); if (dd != NULL) dd->time_frame_duration = temp_float * 1000.; break; case MDC_CONC_BLOCK_SCALE_FACTOR: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tScale factor:\t\t%9.8f\n",temp_float); first_plane->quant_scale = temp_float; break; case MDC_CONC_BLOCK_MINIMUM: /* don't store, we recalculate it */ temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tMinimum:\t\t%9.8f\n",temp_float); break; case MDC_CONC_BLOCK_MAXIMUM: /* don't store, we recalculate it */ temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tMaximum:\t\t%9.8f\n",temp_float); break; case MDC_CONC_BLOCK_DEADTIME_CORRECTION: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tDeadtime correction:\t%9.8f\n",temp_float); break; case MDC_CONC_BLOCK_DECAY_CORRECTION: temp_float = conc_get_float(line, &return_code); if (MDC_INFO) MdcPrntScrn("\tDecay correction:\t%9.8f\n",temp_float); break; case MDC_CONC_BLOCK_PROMPTS: break; /* don't care */ case MDC_CONC_BLOCK_DELAYS: break; /* don't care */ case MDC_CONC_BLOCK_TRUES: break; /* don't care */ case MDC_CONC_BLOCK_PROMPTS_RATE: break; /* don't care */ case MDC_CONC_BLOCK_DELAYS_RATE: break; /* don't care */ case MDC_CONC_BLOCK_SINGLES: break; /* don't care */ case MDC_CONC_BLOCK_END_OF_HEADER: done = MDC_TRUE; break; case MDC_CONC_BLOCK_EOF: done = MDC_TRUE; return("CONC Got inapproprate EOF on reading Concorde/uPET header"); break; case MDC_CONC_BLOCK_UNKNOWN: default: if (num_garbage_lines < MDC_MAX_NUM_GARBAGE_LINES) { MdcPrntWarn("CONC Uninterpretable line: %s",line); } num_garbage_lines++; break; } MdcFree(line); } plane_bytes = first_plane->width*first_plane->height*MdcType2Bytes(first_plane->type); /* where in the raw data file the frame starts */ if ((high_file_pointer >= 0) && (low_file_pointer >= 0)) { #ifdef LONG_8BYTE first_plane->load_location = high_file_pointer*INT_MAX + low_file_pointer; #else /* LONG_4BYTE */ /* since the fseek command uses long, we can't use the value of high_file_pointer in calculating where to fseek. */ if (high_file_pointer != 0) failed_to_read_64_bit = MDC_YES; first_plane->load_location = low_file_pointer; #endif } else { /* DATA_FILE_POINTER header entry not used, calculate it instead */ first_plane->load_location = img*plane_bytes; } /* and read in the info for the first plane, and the rest of the planes */ for(i_plane=0; i_plane < fi->dim[3]; i_plane++, img++) { plane = &fi->image[img]; plane->calibr_fctr = calibration_factor; if (i_plane != 0) { plane->width = first_plane->width; plane->height = first_plane->height; plane->bits = first_plane->bits; plane->type = first_plane->type; plane->quant_scale = first_plane->quant_scale; plane->slice_width = first_plane->slice_width; plane->pixel_xsize = first_plane->pixel_xsize; plane->pixel_ysize = first_plane->pixel_ysize; plane->recon_scale = first_plane->recon_scale; plane->load_location = first_plane->load_location+i_plane*plane_bytes; } } } /* i_frame */ } /* i_gate */ /* complain some more about uninterpretable lines */ if (num_garbage_lines >= MDC_MAX_NUM_GARBAGE_LINES) { MdcPrntWarn("CONC Couldn't process %d header lines", num_garbage_lines); } #ifdef LONG_4BYTE /* warn about trying to access 64 bit locations */ if (failed_to_read_64_bit) { MdcPrntWarn("CONC Read past 2GB (32bit system), file read incomplete"); } #endif return NULL; } const char * MdcLoadCONC(FILEINFO *fi) { const char *msg; msg = MdcLoadHeaderCONC(fi); return(msg); } const char * MdcSavePlaneCONC(FILEINFO *fi, int img) { Int8 saved_norm_over_frames; Uint8 *newbuff, *buff; Int16 pixtype; size_t pixels, bytes; saved_norm_over_frames = MDC_NORM_OVER_FRAMES; if (MDC_QUANTIFY || MDC_CALIBRATE) { /* using global scale factor <=> normalize over ALL images */ MDC_NORM_OVER_FRAMES = MDC_NO; } pixtype = conc_save_type(fi); switch(pixtype) { case BIT16_S: newbuff = MdcGetImgBIT16_S(fi,(unsigned)img); break; case BIT32_S: newbuff = MdcGetImgBIT32_S(fi,(unsigned)img); break; case FLT32: default: newbuff = MdcGetImgFLT32(fi,(unsigned)img); break; } MDC_NORM_OVER_FRAMES = saved_norm_over_frames; if (fi->diff_size == MDC_YES) { buff = MdcGetResizedImage(fi, newbuff, pixtype, (unsigned)img); if (buff == NULL) return("CONC Bad malloc resized image"); MdcFree(newbuff); } else buff = newbuff; if (MDC_FILE_ENDIAN != MDC_HOST_ENDIAN) { MdcMakeImgSwapped(buff,fi,(unsigned)img,fi->mwidth,fi->mheight,pixtype); } pixels = fi->mwidth * fi->mheight; bytes = MdcType2Bytes(pixtype); if (fwrite(buff,bytes,pixels,fi->ofp_raw) != pixels) return("CONC Bad writing of image"); MdcFree(buff); return NULL; } const char *MdcSaveHeaderCONC(FILEINFO *fi) { char raw_filename[MDC_INPUT_STRING_SIZE]; MdcConcFilterTypes filter_type, i_filter_type; IMG_DATA * first_plane; BED_DATA * bd = NULL; GATED_DATA * gd = NULL; int i_bed, i_gate, i_frame, i_plane; float calibration_factor, fstart, fduration, slice_width; Int32 high_file_pointer, low_file_pointer; Uint32 img, fnr; Int16 pixtype; int dimensions, i_dim; size_t write_length; struct tm scan_time; char * pfilename; if (MDC_FILE_STDOUT == MDC_YES) return("CONC Writing to stdout unsupported for this format"); MDC_FILE_ENDIAN = MDC_WRITE_ENDIAN; if (XMDC_GUI == MDC_NO) { MdcDefaultName(fi,MDC_FRMT_CONC,fi->ofname,fi->ifname); } if (MDC_VERBOSE) MdcPrntMesg("Concorde/uPET Writing <%s> ...",fi->ofname); /* check for colored files */ if (fi->map == MDC_MAP_PRESENT) return("CONC Colored files unsupported"); if (MdcKeepFile(fi->ofname)) return("CONC Header file exists!!"); if (fi->dim[7] > 1) return("CONC cannot handle files of this dimensions"); if ((fi->ofp = fopen(fi->ofname, "w")) == NULL) return("CONC Could not open header file for writing"); strncpy(raw_filename,fi->ofname,MDC_INPUT_STRING_SIZE-5); pfilename = strstr(raw_filename, ".img.hdr"); if (pfilename != NULL) strcpy(pfilename+4, "\0"); else strcat(raw_filename,".dat"); if (MdcKeepFile(raw_filename)) return("CONC Image file exists!!"); if ((fi->ofp_raw = fopen(raw_filename, "wb")) == NULL) return("CONC Could not open data file for writing"); fprintf(fi->ofp, "#\n# Header file for data file %s\n", raw_filename); fprintf(fi->ofp, "#\twith %d frames\n",fi->dim[4]*fi->dim[5]); fprintf(fi->ofp, "#\n# Concorde/µPET image file - %s %s\n#\n" ,XMEDCON_PRGR, XMEDCON_VERSION); fprintf(fi->ofp, "#\n%s %5.3f\n", MdcConcHdrValueNames[MDC_CONC_HDR_VERSION] , MDC_CONC_SUPPORTED_VERSION); fprintf(fi->ofp, "#\n%s %s\n", MdcConcHdrValueNames[MDC_CONC_HDR_INSTITUTION] , fi->institution); fprintf(fi->ofp, "#\n%s %s\n", MdcConcHdrValueNames[MDC_CONC_HDR_STUDY] , fi->study_id); fprintf(fi->ofp, "#\n%s %s\n", MdcConcHdrValueNames[MDC_CONC_HDR_FILE_NAME] , raw_filename); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_FILE_TYPE] , MDC_CONC_FILE_IMAGE); switch(fi->acquisition_type) { case MDC_ACQUISITION_TOMO: case MDC_ACQUISITION_STATIC: fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ACQUISITION_MODE] , MDC_CONC_ACQ_EMISSION); break; case MDC_ACQUISITION_DYNAMIC: fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ACQUISITION_MODE] , MDC_CONC_ACQ_DYNAMIC); break; case MDC_ACQUISITION_GATED: fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ACQUISITION_MODE] , MDC_CONC_ACQ_GATED); break; case MDC_ACQUISITION_GSPECT: case MDC_ACQUISITION_UNKNOWN: default: fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ACQUISITION_MODE] , MDC_CONC_ACQ_UNKNOWN); break; } fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_TOTAL_FRAMES] , fi->dim[4]*fi->dim[5]); fprintf(fi->ofp, "#\n%s %s\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ISOTOPE] , fi->isotope_code); fprintf(fi->ofp, "#\n%s %e\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ISOTOPE_HALF_LIFE] , fi->isotope_halflife); fprintf(fi->ofp, "# Note: isotope branching fraction is included in the calibration fraction\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_ISOTOPE_BRANCHING_FRACTION] , 1.0); slice_width = fi->pixdim[3]; #ifdef MDC_USE_SLICE_SPACING if (fi->number > 1) slice_width = fi->image[0].slice_spacing; #endif fprintf(fi->ofp, "#\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_AXIAL_CRYSTAL_PITCH] , 2.0*slice_width/10.0); pixtype = conc_save_type(fi); switch (pixtype) { case BIT8_S: fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_SBYTE); break; case BIT16_S: if (MDC_FILE_ENDIAN == MDC_LITTLE_ENDIAN) { fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_SSHORT_LE); }else{ fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_SSHORT_BE); } break; case BIT32_S: if (MDC_FILE_ENDIAN == MDC_LITTLE_ENDIAN) { fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_SINT_LE); }else{ fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_SINT_BE); } break; case FLT32 : default: if (MDC_FILE_ENDIAN == MDC_LITTLE_ENDIAN) { fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_FLOAT_LE); }else{ fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_TYPE] , MDC_CONC_DATA_FLOAT_BE); } } fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_DATA_ORDER] , MDC_CONC_ORDER_SINOGRAM); dimensions = 0; for (i_dim=1;i_dim<=6;i_dim++) dimensions += (fi->dim[i_dim] > 1) ? 1 : 0; fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_NUMBER_OF_DIMENSIONS] , 3); /* dimension is always 3... regardless of dynamic or gated */ fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_X_DIMENSION] , fi->dim[1]); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_Y_DIMENSION] , fi->dim[2]); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_Z_DIMENSION] , fi->dim[3]); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_W_DIMENSION] , 1); /* not sure what this is used for... */ filter_type = 0; for (i_filter_type=0;i_filter_typefilter_type) == 0) filter_type = i_filter_type; fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_X_FILTER] , filter_type); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_Y_FILTER] , MDC_CONC_FILTER_NONE); fprintf(fi->ofp, "#\n%s %d\n", MdcConcHdrValueNames[MDC_CONC_HDR_Z_FILTER] , MDC_CONC_FILTER_NONE); fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_RECON_ALGORITHM] , MDC_CONC_RECON_UNKNOWN); fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_DECAY_CORRECTION_APPLIED] , fi->decay_corrected); fprintf(fi->ofp, "#\n%s %g\n", MdcConcHdrValueNames[MDC_CONC_HDR_PIXEL_SIZE] , fi->pixdim[1]/10.0); /* eNlf: don't use, scales always combined internally calibration_factor = fi->image[0].calibr_fctr; */ calibration_factor = 1.0; fprintf(fi->ofp, "#\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_CALIBRATION_FACTOR] , calibration_factor); if ((fi->study_date_month != 0) && (fi->study_date_year != 0)) { scan_time.tm_sec = fi->study_time_second; scan_time.tm_min = fi->study_time_minute; scan_time.tm_hour = fi->study_time_hour; scan_time.tm_mday = fi->study_date_day; scan_time.tm_mon = fi->study_date_month-1; scan_time.tm_year = fi->study_date_year-1900; scan_time.tm_isdst = -1; /* "-1" is suppose to let the system figure it out */ if (mktime(&scan_time) != -1) { /* make sure the time is proper */ fprintf(fi->ofp, "#\n%s %s" , MdcConcHdrValueNames[MDC_CONC_HDR_SCAN_TIME] , asctime(&scan_time)); } } fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_DOSE_UNITS] , MDC_CONC_DOSE_UNITS_MEGA_BEQUERELS); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_INJECTED_DOSE] , fi->injected_dose); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_INJECTION_DECAY_CORRECTION] , 1.0); fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_GATE_INPUTS] , fi->gatednr); for (i_gate=0; i_gate < fi->gatednr; i_gate++) { gd = &(fi->gdata[i_gate]); fprintf(fi->ofp, "#\n%s %d %1.0f %g %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_GATE_BINS] , i_gate , gd->nr_projections , gd->window_low/1000. , gd->window_high/1000.); } fprintf(fi->ofp, "#\n%s %d\n" , MdcConcHdrValueNames[MDC_CONC_HDR_SUBJECT_WEIGHT_UNITS] , MDC_CONC_WEIGHT_UNITS_KILOGRAMS); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcHdrValueNames[MDC_CONC_HDR_SUBJECT_WEIGHT] , fi->patient_weight); fprintf(fi->ofp, "#\n%s\n", MdcConcHdrValueNames[MDC_CONC_HDR_END_OF_HEADER]); /* write the data and the headers for the frames */ fprintf(fi->ofp, "#\n#\n#\n#\n"); img = 0; high_file_pointer = 0; low_file_pointer = 0; for (i_bed = 0; i_bed < fi->dim[6]; i_bed++) { if (fi->bednr > 0) bd = &fi->beddata[i_bed]; for (i_gate = 0; i_gate < fi->dim[5]; i_gate++) { for (i_frame= 0; i_frame < fi->dim[4]; i_frame++) { first_plane = &fi->image[img]; fnr = first_plane->frame_number; if ((fi->dynnr > 0) && (fnr > 0)) { fstart = fi->dyndata[fnr-1].time_frame_start/1000.; fduration = fi->dyndata[fnr-1].time_frame_duration/1000.; }else{ fstart = 0.; fduration = 0.; } fprintf(fi->ofp, "#\n%s %d\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_FRAME], i_frame+i_gate*fi->dim[4]); /* Concorde's ASIPro program requires the event_type entry on gated data for some reason... */ fprintf(fi->ofp, "#\n%s %d\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_EVENT_TYPE], MDC_CONC_EVENT_UNKNOWN); fprintf(fi->ofp, "#\n%s %d\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_GATE], i_gate); fprintf(fi->ofp, "#\n%s %d\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_BED], i_bed); if (bd != NULL) { fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_BED_OFFSET] , bd->hoffset/10.); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_VERTICAL_BED_OFFSET] , bd->voffset/10.); } fprintf(fi->ofp, "#\n#\tData file offset to start of data," \ " two 32 bit signed ints\n"); fprintf(fi->ofp, "%s %d %d\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_DATA_FILE_POINTER] , high_file_pointer, low_file_pointer); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_FRAME_START] , fstart); fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_FRAME_DURATION] , fduration); for (i_plane = 0; i_plane < fi->dim[3]; i_plane++, img++) { write_length = fi->mwidth * fi->mheight*MdcType2Bytes(pixtype); /* update the high and low file pointers */ if ((INT_MAX-write_length) < low_file_pointer) { high_file_pointer += 1; low_file_pointer = write_length - (INT_MAX-low_file_pointer); } else low_file_pointer += write_length; } if (first_plane->rescaled) { fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_SCALE_FACTOR] , first_plane->rescaled_fctr); } else { /* eNlf: must use combined scale factor fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_SCALE_FACTOR] , first_plane->quant_scale); */ fprintf(fi->ofp, "#\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_SCALE_FACTOR] , first_plane->rescale_slope); } /* Concorde's ASIPro program requires the dead time correction entry on gated data.... */ fprintf(fi->ofp, "#\n# Not 1.0, Unknown\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_DEADTIME_CORRECTION] , 1.0); /* Concorde's ASIPro program requires the decay correction entry on gated data.... */ { float num_half_lifes = 0.; if (fi->isotope_halflife > 0.) { num_half_lifes = (fstart+fduration/2.0)/fi->isotope_halflife; } fprintf(fi->ofp, "#\n# Check decay_correction_applied to know if already applied\n%s %g\n" , MdcConcBlockValueNames[MDC_CONC_BLOCK_DECAY_CORRECTION] , pow(0.5,(double)num_half_lifes)); } fprintf(fi->ofp, "#\n%s\n", MdcConcBlockValueNames[MDC_CONC_BLOCK_END_OF_HEADER]); } } } return NULL; } const char *MdcSaveCONC(FILEINFO *fi) { const char * return_string; int img=0; int i_bed, i_gate, i_frame, i_plane; return_string = MdcSaveHeaderCONC(fi); if (return_string != NULL) return(return_string); for (i_bed = 0; i_bed < fi->dim[6]; i_bed++) { for (i_gate = 0; i_gate < fi->dim[5]; i_gate++) { for (i_frame= 0; i_frame < fi->dim[4]; i_frame++) { for (i_plane = 0; i_plane < fi->dim[3]; i_plane++, img++) { return_string = MdcSavePlaneCONC(fi, img); if (return_string != NULL) return(return_string); } } } } MdcCheckQuantitation(fi); return(NULL); } int MdcCheckCONC(FILEINFO *fi) { char header_str[17]; int FORMAT=MDC_FRMT_NONE; if (fscanf(fi->ifp, "%16s", header_str)==0) return(MDC_BAD_READ); if (strcmp(header_str, "#") == 0) { if (fscanf(fi->ifp, "%16s", header_str)==0) return(MDC_BAD_READ); if (strcmp(header_str, "#") == 0) { if (fscanf(fi->ifp, "%16s", header_str)==0) return(MDC_BAD_READ); if (strcmp(header_str, "Header") == 0) { if (fscanf(fi->ifp, "%16s", header_str)==0) return(MDC_BAD_READ); if (strcmp(header_str, "file") == 0) return (MDC_FRMT_CONC); } } } return(FORMAT); } const char *MdcReadCONC(FILEINFO *fi) { const char * return_string; int i_gate, i_frame, i_plane, img=0; int total_images=0; if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_BEGIN,0.,"Reading Concorde/uPET:"); /* read in the header */ return_string = MdcLoadHeaderCONC(fi); if (return_string != NULL) return(return_string); total_images = fi->dim[4]*fi->dim[3]; /* make sure all planes are loaded */ for (i_gate = 0; i_gate < fi->dim[5]; i_gate++) { for (i_frame = 0; i_frame < fi->dim[4]; i_frame++) { if (MDC_PROGRESS && (total_images > 100)) { MdcProgress(MDC_PROGRESS_INCR,1./(float)(fi->dim[4]*fi->dim[5]),NULL); } /* and read in the data for the first plane, and the rest of the planes */ for (i_plane=0; i_plane < fi->dim[3]; i_plane++, img++) { if (MDC_PROGRESS && (total_images <= 100)) { MdcProgress(MDC_PROGRESS_INCR,1./(float)fi->dim[3],NULL); } return_string = MdcLoadPlaneCONC(fi, img); if (return_string != NULL) return(return_string); } } } return(NULL); } const char *MdcWriteCONC(FILEINFO *fi) { const char * return_string; int img=0; int i_bed, i_gate, i_frame, i_plane; int total_images=0; if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_BEGIN,0.,"Writing Concorde/uPET:"); total_images = fi->dim[4]*fi->dim[3]; return_string = MdcSaveHeaderCONC(fi); if (return_string != NULL) return(return_string); for (i_bed = 0; i_bed < fi->dim[6]; i_bed++) { for (i_gate = 0; i_gate < fi->dim[5]; i_gate++) { for (i_frame= 0; i_frame < fi->dim[4]; i_frame++) { if (MDC_PROGRESS && (total_images > 100)) { MdcProgress(MDC_PROGRESS_INCR,1./(float)fi->dim[4],NULL); } for (i_plane = 0; i_plane < fi->dim[3]; i_plane++, img++) { if (MDC_PROGRESS && (total_images <= 100)) { MdcProgress(MDC_PROGRESS_INCR,1./(float)fi->dim[3],NULL); } return_string = MdcSavePlaneCONC(fi, img); if (return_string != NULL) return(return_string); } } } } MdcCheckQuantitation(fi); return(NULL); }