From 4472c115d5a14faa284334517645a911a101847b Mon Sep 17 00:00:00 2001 From: Jonas Teuwen Date: Sat, 22 Apr 2023 13:40:58 +0200 Subject: [PATCH] Implement read_region and isyntax-to-tiff writer --- CMakeLists.txt | 73 +- src/isyntax/isyntax.c | 5080 +++++++++++++++++----------------- src/isyntax/isyntax_reader.c | 3 +- src/isyntax_example.c | 2 +- src/isyntax_to_tiff.c | 422 +++ src/libisyntax.c | 90 + src/libisyntax.h | 2 + 7 files changed, 3126 insertions(+), 2546 deletions(-) create mode 100644 src/isyntax_to_tiff.c diff --git a/CMakeLists.txt b/CMakeLists.txt index 43f9ca3..afb6008 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -4,15 +4,53 @@ project(libisyntax) set(CMAKE_C_STANDARD 17) if(NOT CMAKE_BUILD_TYPE) - set(CMAKE_BUILD_TYPE Release) + set(CMAKE_BUILD_TYPE Debug) +endif() + +if(CMAKE_SYSTEM_PROCESSOR MATCHES "^(aarch64.*|AARCH64.*|arm64.*|ARM64.*)") + set(IS_ARM64 TRUE) + set(ARM_ARCH "arm64") +elseif(CMAKE_SYSTEM_PROCESSOR MATCHES "^(arm.*|ARM.*)") + set(IS_ARM TRUE) + set(ARM_ARCH "arm") endif() -message(STATUS "CMAKE_BUILD_TYPE = ${CMAKE_BUILD_TYPE}") -if(CMAKE_SYSTEM_PROCESSOR MATCHES "^(arm|aarch64)") - message(STATUS "Enabling ARM NEON support") - add_compile_options(-mfpu=neon) +if(APPLE AND (IS_ARM OR IS_ARM64)) + set(IS_APPLE_SILICON TRUE) endif() +if(IS_APPLE_SILICON) + message(STATUS "Detected Apple silicon ${ARM_ARCH} architecture") + set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -arch ${ARM_ARCH}") + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -arch ${ARM_ARCH}") + # Enable NEON for ARM-based Apple Silicon + add_compile_options(-march=armv8.2-a+fp16+simd) +endif() + + +option(ENABLE_MEMORY_CHECK "Enable memory check with sanitizers (requires DEBUG mode)" OFF) + +if(ENABLE_MEMORY_CHECK) + if(NOT CMAKE_BUILD_TYPE MATCHES Debug) + message(FATAL_ERROR "ENABLE_MEMORY_CHECK is only allowed in Debug build mode") + endif() + + if(APPLE) + set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -fsanitize=address -O1 -fno-omit-frame-pointer -g") + set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} -fsanitize=address -fno-omit-frame-pointer") + set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -fsanitize=address -fno-omit-frame-pointer") + set(CMAKE_EXE_LINKER_FLAGS_DEBUG "${CMAKE_EXE_LINKER_FLAGS_DEBUG} -fsanitize=address") + elseif(UNIX AND NOT APPLE) # Linux + set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -O1 -g") + set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} -Og -g") + set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -Og -g") + endif() +endif() + +message(STATUS "CMAKE_BUILD_TYPE = ${CMAKE_BUILD_TYPE}") + + + set(CMAKE_RUNTIME_OUTPUT_DIRECTORY "${CMAKE_SOURCE_DIR}") include_directories("${CMAKE_SOURCE_DIR}/src") @@ -21,6 +59,15 @@ include_directories("${CMAKE_SOURCE_DIR}/src/utils") include_directories("${CMAKE_SOURCE_DIR}/src/isyntax") include_directories("${CMAKE_SOURCE_DIR}/src/third_party") + +# Add threads library for linux +find_package(Threads REQUIRED) + +if(UNIX AND NOT APPLE) + set(LIBRT rt) +endif() + + set(LIBISYNTAX_COMMON_SOURCE_FILES src/libisyntax.c src/isyntax/isyntax.c @@ -49,6 +96,22 @@ add_executable(isyntax_example ${LIBISYNTAX_COMMON_SOURCE_FILES} ) +target_link_libraries(isyntax_example PRIVATE Threads::Threads ${LIBRT}) + +# Find LibTIFF library +find_package(TIFF REQUIRED) +if (NOT TIFF_FOUND) + message(WARNING "LibTIFF not found") + message(WARNING "Will not compile `isyntax-to-tiff` utility") +else() + include_directories(${TIFF_INCLUDE_DIR}) + add_executable(isyntax-to-tiff + src/isyntax_to_tiff.c + ${LIBISYNTAX_COMMON_SOURCE_FILES} ${TIFF_LIBRARIES} src/isyntax_to_tiff.c) + target_link_libraries(isyntax-to-tiff PRIVATE ${TIFF_LIBRARIES} Threads::Threads ${LIBRT}) +endif() + + if (WIN32) target_link_libraries(libisyntax winmm) target_link_libraries(isyntax_example winmm) diff --git a/src/isyntax/isyntax.c b/src/isyntax/isyntax.c index 765b81b..2839a83 100644 --- a/src/isyntax/isyntax.c +++ b/src/isyntax/isyntax.c @@ -95,238 +95,238 @@ * See README for more details. */ static const unsigned char base64_table[65] = - "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"; + "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"; unsigned char * base64_decode(const unsigned char *src, size_t len, size_t *out_len) { - unsigned char dtable[256], *out, *pos, block[4], tmp; - size_t i, count, olen; - int pad = 0; - - memset(dtable, 0x80, 256); - for (i = 0; i < sizeof(base64_table) - 1; i++) - dtable[base64_table[i]] = (unsigned char) i; - dtable['='] = 0; - - count = 0; - for (i = 0; i < len; i++) { - if (dtable[src[i]] != 0x80) - count++; - } + unsigned char dtable[256], *out, *pos, block[4], tmp; + size_t i, count, olen; + int pad = 0; + + memset(dtable, 0x80, 256); + for (i = 0; i < sizeof(base64_table) - 1; i++) + dtable[base64_table[i]] = (unsigned char) i; + dtable['='] = 0; + + count = 0; + for (i = 0; i < len; i++) { + if (dtable[src[i]] != 0x80) + count++; + } - if (count == 0 || count % 4) - return NULL; - - olen = count / 4 * 3; - pos = out = malloc(olen); - if (out == NULL) - return NULL; - - count = 0; - for (i = 0; i < len; i++) { - tmp = dtable[src[i]]; - if (tmp == 0x80) - continue; - - if (src[i] == '=') - pad++; - block[count] = tmp; - count++; - if (count == 4) { - *pos++ = (block[0] << 2) | (block[1] >> 4); - *pos++ = (block[1] << 4) | (block[2] >> 2); - *pos++ = (block[2] << 6) | block[3]; - count = 0; - if (pad) { - if (pad == 1) - pos--; - else if (pad == 2) - pos -= 2; - else { - /* Invalid padding */ - free(out); - return NULL; - } - break; - } - } - } + if (count == 0 || count % 4) + return NULL; + + olen = count / 4 * 3; + pos = out = malloc(olen); + if (out == NULL) + return NULL; + + count = 0; + for (i = 0; i < len; i++) { + tmp = dtable[src[i]]; + if (tmp == 0x80) + continue; + + if (src[i] == '=') + pad++; + block[count] = tmp; + count++; + if (count == 4) { + *pos++ = (block[0] << 2) | (block[1] >> 4); + *pos++ = (block[1] << 4) | (block[2] >> 2); + *pos++ = (block[2] << 6) | block[3]; + count = 0; + if (pad) { + if (pad == 1) + pos--; + else if (pad == 2) + pos -= 2; + else { + /* Invalid padding */ + free(out); + return NULL; + } + break; + } + } + } - *out_len = pos - out; - return out; + *out_len = pos - out; + return out; } // end of base64 decoder. // similar to atoi(), but also returning the string position so we can chain calls one after another. static const char* atoi_and_advance(const char* str, i32* dest) { - i32 num = 0; - bool neg = false; - while (isspace(*str)) ++str; - if (*str == '-') { - neg = true; - ++str; - } - while (isdigit(*str)) { - num = 10*num + (*str - '0'); - ++str; - } - if (neg) num = -num; - *dest = num; - return str; + i32 num = 0; + bool neg = false; + while (isspace(*str)) ++str; + if (*str == '-') { + neg = true; + ++str; + } + while (isdigit(*str)) { + num = 10*num + (*str - '0'); + ++str; + } + if (neg) num = -num; + *dest = num; + return str; } static void parse_three_integers(const char* str, i32* first, i32* second, i32* third) { - str = atoi_and_advance(str, first); - str = atoi_and_advance(str, second); - atoi_and_advance(str, third); + str = atoi_and_advance(str, first); + str = atoi_and_advance(str, second); + atoi_and_advance(str, third); } static i32 parse_up_to_five_integers(const char* str, i32* array_of_five_integers) { - for (i32 i = 0; i < 5; ++i) { - if (*str == '\0') { - for (i32 j = i; i < 5; ++j) { - array_of_five_integers[j] = 0; - return i; // return number of valid integers - } - } else { - str = atoi_and_advance(str, array_of_five_integers + i); - } - } - return 5; + for (i32 i = 0; i < 5; ++i) { + if (*str == '\0') { + for (i32 j = i; i < 5; ++j) { + array_of_five_integers[j] = 0; + return i; // return number of valid integers + } + } else { + str = atoi_and_advance(str, array_of_five_integers + i); + } + } + return 5; } static void isyntax_parse_ufsimport_child_node(isyntax_t* isyntax, u32 group, u32 element, char* value, u64 value_len) { - switch(group) { - default: { - // unknown group - console_print_verbose("Unknown group 0x%04x\n", group); - } break; - case 0x0008: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x002A: /*DICOM_ACQUISITION_DATETIME*/ {} break; // "20210101103030.000000" - case 0x0070: /*DICOM_MANUFACTURER*/ {} break; // "PHILIPS" - case 0x1090: /*DICOM_MANUFACTURERS_MODEL_NAME*/ {} break; // "UFS Scanner" - } - }; break; - case 0x0018: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x1000: /*DICOM_DEVICE_SERIAL_NUMBER*/ {} break; // "FMT<4-digit number>" - case 0x1020: /*DICOM_SOFTWARE_VERSIONS*/ {} break; // "" "" - case 0x1200: /*DICOM_DATE_OF_LAST_CALIBRATION*/ {} break; // "20210101" - case 0x1201: /*DICOM_TIME_OF_LAST_CALIBRATION*/ {} break; // "100730" - } - } break; - case 0x101D: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x1007: /*PIIM_DP_SCANNER_RACK_NUMBER*/ {} break; // "[1..15]" - case 0x1008: /*PIIM_DP_SCANNER_SLOT_NUMBER*/ {} break; // "[1..15]" - case 0x1009: /*PIIM_DP_SCANNER_OPERATOR_ID*/ {} break; // "" - case 0x100A: /*PIIM_DP_SCANNER_CALIBRATION_STATUS*/ {} break; // "OK" or "NOT OK" - - } - } break; - case 0x301D: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x1001: /*PIM_DP_UFS_INTERFACE_VERSION*/ { - // Value will likely be "5.0" for v1 iSyntax files, "100.5" for v2 iSyntax files - isyntax->data_model_major_version = atoi(value); - } break; - case 0x1002: /*PIM_DP_UFS_BARCODE*/ {} break; // "" - case 0x1003: /*PIM_DP_SCANNED_IMAGES*/ {} break; - case 0x1010: /*PIM_DP_SCANNER_RACK_PRIORITY*/ {} break; // "" - } - } break; - } + switch(group) { + default: { + // unknown group + console_print_verbose("Unknown group 0x%04x\n", group); + } break; + case 0x0008: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x002A: /*DICOM_ACQUISITION_DATETIME*/ {} break; // "20210101103030.000000" + case 0x0070: /*DICOM_MANUFACTURER*/ {} break; // "PHILIPS" + case 0x1090: /*DICOM_MANUFACTURERS_MODEL_NAME*/ {} break; // "UFS Scanner" + } + }; break; + case 0x0018: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x1000: /*DICOM_DEVICE_SERIAL_NUMBER*/ {} break; // "FMT<4-digit number>" + case 0x1020: /*DICOM_SOFTWARE_VERSIONS*/ {} break; // "" "" + case 0x1200: /*DICOM_DATE_OF_LAST_CALIBRATION*/ {} break; // "20210101" + case 0x1201: /*DICOM_TIME_OF_LAST_CALIBRATION*/ {} break; // "100730" + } + } break; + case 0x101D: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x1007: /*PIIM_DP_SCANNER_RACK_NUMBER*/ {} break; // "[1..15]" + case 0x1008: /*PIIM_DP_SCANNER_SLOT_NUMBER*/ {} break; // "[1..15]" + case 0x1009: /*PIIM_DP_SCANNER_OPERATOR_ID*/ {} break; // "" + case 0x100A: /*PIIM_DP_SCANNER_CALIBRATION_STATUS*/ {} break; // "OK" or "NOT OK" + + } + } break; + case 0x301D: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x1001: /*PIM_DP_UFS_INTERFACE_VERSION*/ { + // Value will likely be "5.0" for v1 iSyntax files, "100.5" for v2 iSyntax files + isyntax->data_model_major_version = atoi(value); + } break; + case 0x1002: /*PIM_DP_UFS_BARCODE*/ {} break; // "" + case 0x1003: /*PIM_DP_SCANNED_IMAGES*/ {} break; + case 0x1010: /*PIM_DP_SCANNER_RACK_PRIORITY*/ {} break; // "" + } + } break; + } } static bool isyntax_parse_scannedimage_child_node(isyntax_t* isyntax, u32 group, u32 element, char* value, u64 value_len) { - // Parse metadata belong to one of the images in the file (either a WSI, LABELIMAGE or MACROIMAGE) + // Parse metadata belong to one of the images in the file (either a WSI, LABELIMAGE or MACROIMAGE) - isyntax_image_t* image = isyntax->parser.current_image; - if (!image) { - image = isyntax->parser.current_image = &isyntax->images[0]; - } + isyntax_image_t* image = isyntax->parser.current_image; + if (!image) { + image = isyntax->parser.current_image = &isyntax->images[0]; + } - bool success = true; - - switch(group) { - default: { - // unknown group - console_print_verbose("Unknown group 0x%04x\n", group); - } break; - case 0x0008: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x2111: /*DICOM_DERIVATION_DESCRIPTION*/ { // "PHILIPS UFS V%s | Quality=%d | DWT=%d | Compressor=%d" - - } break; - } - }; break; - case 0x0028: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x0002: /*DICOM_SAMPLES_PER_PIXEL*/ {} break; - case 0x0100: /*DICOM_BITS_ALLOCATED*/ {} break; - case 0x0101: /*DICOM_BITS_STORED*/ {} break; - case 0x0102: /*DICOM_HIGH_BIT*/ {} break; - case 0x0103: /*DICOM_PIXEL_REPRESENTATION*/ {} break; - case 0x2000: /*DICOM_ICCPROFILE*/ {} break; - case 0x2110: /*DICOM_LOSSY_IMAGE_COMPRESSION*/ {} break; - case 0x2112: /*DICOM_LOSSY_IMAGE_COMPRESSION_RATIO*/ {} break; - case 0x2114: /*DICOM_LOSSY_IMAGE_COMPRESSION_METHOD*/ {} break; // "PHILIPS_DP_1_0" - } - } break; - case 0x301D: { - switch(element) { - default: { - console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); - } break; - case 0x1004: /*PIM_DP_IMAGE_TYPE*/ { // "MACROIMAGE" or "LABELIMAGE" or "WSI" - if ((strcmp(value, "MACROIMAGE") == 0)) { - isyntax->macro_image_index = isyntax->parser.running_image_index; - isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_MACROIMAGE; - image->image_type = ISYNTAX_IMAGE_TYPE_MACROIMAGE; - } else if ((strcmp(value, "LABELIMAGE") == 0)) { - isyntax->label_image_index = isyntax->parser.running_image_index; - isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_LABELIMAGE; - image->image_type = ISYNTAX_IMAGE_TYPE_LABELIMAGE; - } else if ((strcmp(value, "WSI") == 0)) { - isyntax->wsi_image_index = isyntax->parser.running_image_index; - isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_WSI; - image->image_type = ISYNTAX_IMAGE_TYPE_WSI; - } - } break; - case 0x1005: { /*PIM_DP_IMAGE_DATA*/ - size_t decoded_capacity = value_len; - size_t decoded_len = 0; - i32 last_char = value[value_len-1]; - if (last_char == '/') { - value_len--; // The last character may cause the base64 decoding to fail if invalid - } - u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); - if (decoded) { - i32 channels_in_file = 0; + bool success = true; + + switch(group) { + default: { + // unknown group + console_print_verbose("Unknown group 0x%04x\n", group); + } break; + case 0x0008: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x2111: /*DICOM_DERIVATION_DESCRIPTION*/ { // "PHILIPS UFS V%s | Quality=%d | DWT=%d | Compressor=%d" + + } break; + } + }; break; + case 0x0028: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x0002: /*DICOM_SAMPLES_PER_PIXEL*/ {} break; + case 0x0100: /*DICOM_BITS_ALLOCATED*/ {} break; + case 0x0101: /*DICOM_BITS_STORED*/ {} break; + case 0x0102: /*DICOM_HIGH_BIT*/ {} break; + case 0x0103: /*DICOM_PIXEL_REPRESENTATION*/ {} break; + case 0x2000: /*DICOM_ICCPROFILE*/ {} break; + case 0x2110: /*DICOM_LOSSY_IMAGE_COMPRESSION*/ {} break; + case 0x2112: /*DICOM_LOSSY_IMAGE_COMPRESSION_RATIO*/ {} break; + case 0x2114: /*DICOM_LOSSY_IMAGE_COMPRESSION_METHOD*/ {} break; // "PHILIPS_DP_1_0" + } + } break; + case 0x301D: { + switch(element) { + default: { + console_print_verbose("Unknown element (0x%04x, 0x%04x)\n", group, element); + } break; + case 0x1004: /*PIM_DP_IMAGE_TYPE*/ { // "MACROIMAGE" or "LABELIMAGE" or "WSI" + if ((strcmp(value, "MACROIMAGE") == 0)) { + isyntax->macro_image_index = isyntax->parser.running_image_index; + isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_MACROIMAGE; + image->image_type = ISYNTAX_IMAGE_TYPE_MACROIMAGE; + } else if ((strcmp(value, "LABELIMAGE") == 0)) { + isyntax->label_image_index = isyntax->parser.running_image_index; + isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_LABELIMAGE; + image->image_type = ISYNTAX_IMAGE_TYPE_LABELIMAGE; + } else if ((strcmp(value, "WSI") == 0)) { + isyntax->wsi_image_index = isyntax->parser.running_image_index; + isyntax->parser.current_image_type = ISYNTAX_IMAGE_TYPE_WSI; + image->image_type = ISYNTAX_IMAGE_TYPE_WSI; + } + } break; + case 0x1005: { /*PIM_DP_IMAGE_DATA*/ + size_t decoded_capacity = value_len; + size_t decoded_len = 0; + i32 last_char = value[value_len-1]; + if (last_char == '/') { + value_len--; // The last character may cause the base64 decoding to fail if invalid + } + u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); + if (decoded) { + i32 channels_in_file = 0; #ifdef ISYNTAX_JPEG_DECODER_USE_LIBJPEG - // TODO: Why does this crash? + // TODO: Why does this crash? // Apparently, there is a bug in the libjpeg-turbo implementation of jsimd_can_h2v2_fancy_upsample() when using SIMD. // jsimd_h2v2_fancy_upsample_avx2 writes memory out of bounds. // This causes the program to crash eventually when trying to free memory in free_pool(). @@ -336,880 +336,883 @@ static bool isyntax_parse_scannedimage_child_node(isyntax_t* isyntax, u32 group, image->pixels = jpeg_decode_image(decoded, decoded_len, &image->width, &image->height, &channels_in_file); #else - // stb_image.h - image->pixels = stbi_load_from_memory(decoded, decoded_len, &image->width, &image->height, &channels_in_file, 4); + // stb_image.h + image->pixels = stbi_load_from_memory(decoded, decoded_len, &image->width, &image->height, &channels_in_file, 4); #endif - free(decoded); - } - } break; - case 0x1013: /*DP_COLOR_MANAGEMENT*/ {} break; - case 0x1014: /*DP_IMAGE_POST_PROCESSING*/ {} break; - case 0x1015: /*DP_SHARPNESS_GAIN_RGB24*/ {} break; - case 0x1016: /*DP_CLAHE_CLIP_LIMIT_Y16*/ {} break; - case 0x1017: /*DP_CLAHE_NR_BINS_Y16*/ {} break; - case 0x1018: /*DP_CLAHE_CONTEXT_DIMENSION_Y16*/ {} break; - case 0x1019: /*DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR*/ {} break; - case 0x101A: /*DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL*/ {} break; - case 0x101B: /*DP_WAVELET_QUANTIZER*/ {} break; - case 0x101C: /*DP_WAVELET_DEADZONE*/ {} break; - case 0x1025: /*UFS_IMAGE_OPP_EXTREME_VERTEX*/ {} break; // data model >= 100 - case 0x2000: /*UFS_IMAGE_GENERAL_HEADERS*/ {} break; - case 0x2001: /*UFS_IMAGE_NUMBER_OF_BLOCKS*/ { - // NOTE: the actual number of stored codeblocks may be lower than this number. - image->number_of_blocks = atoi(value); - } break; - case 0x2002: /*UFS_IMAGE_DIMENSIONS_OVER_BLOCK*/ {} break; - case 0x2003: /*UFS_IMAGE_DIMENSIONS*/ {} break; - case 0x2004: /*UFS_IMAGE_DIMENSION_NAME*/ {} break; - case 0x2005: /*UFS_IMAGE_DIMENSION_TYPE*/ {} break; - case 0x2006: /*UFS_IMAGE_DIMENSION_UNIT*/ {} break; - case 0x2007: /*UFS_IMAGE_DIMENSION_SCALE_FACTOR*/ { - float mpp = atof(value); - if (isyntax->parser.dimension_index == 0 /*x*/) { - isyntax->mpp_x = mpp; - isyntax->is_mpp_known = true; - } else if (isyntax->parser.dimension_index == 1 /*y*/) { - isyntax->mpp_y = mpp; - isyntax->is_mpp_known = true; - } - } break; - case 0x2008: /*UFS_IMAGE_DIMENSION_DISCRETE_VALUES_STRING*/ {} break; - case 0x2009: /*UFS_IMAGE_BLOCK_HEADER_TEMPLATES*/ {} break; - case 0x200A: /*UFS_IMAGE_DIMENSION_RANGES*/ {} break; - case 0x200B: /*UFS_IMAGE_DIMENSION_RANGE*/ { - isyntax_image_dimension_range_t range = {0}; - parse_three_integers(value, &range.start, &range.step, &range.end); - i32 step_nonzero = (range.step != 0) ? range.step : 1; - range.numsteps = ((range.end + range.step) - range.start) / step_nonzero; - if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate) { - isyntax_block_header_template_t* template = isyntax->block_header_templates + isyntax->parser.block_header_template_index; - switch(isyntax->parser.dimension_index) { - default: break; - case 0: template->block_width = range.numsteps; break; - case 1: template->block_height = range.numsteps; break; - case 2: template->color_component = range.start; break; - case 3: template->scale = range.start; break; - case 4: template->waveletcoeff = (range.start == 0) ? 1 : 3; break; - } - DUMMY_STATEMENT; - } else if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageGeneralHeader) { - switch(isyntax->parser.dimension_index) { - default: break; - case 0: { - image->offset_x = range.start; - image->width = range.numsteps; - } break; - case 1: { - image->offset_y = range.start; - image->height = range.numsteps; - } break; - case 2: break; // always 3 color channels ("Y" "Co" "Cg"), no need to check - case 3: { - image->level_count = range.numsteps; - image->max_scale = range.numsteps - 1; - } break; - case 4: break; // always 4 wavelet coefficients ("LL" "LH" "HL" "HH"), no need to check - } - DUMMY_STATEMENT; - } else if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate) { - // The dimension range information within a cluster header template (specifically, the range.start) - // seems to provide the 'base values' for the dimension offsets for each of the codeblocks in the cluster template. - // To get the correct coordinate info for each codeblock, we need to add these 'base values' - // to the codeblock-specific coordinate values from the UFS_IMAGE_BLOCK_COORDINATE tags. - isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; - switch(isyntax->parser.dimension_index) { - default: break; - case 0: template->base_x = range.start; break; - case 1: template->base_y = range.start; break; - case 2: template->base_color_component = range.start; break; - case 3: template->base_scale = range.start; break; - case 4: template->base_waveletcoeff = range.start; break; - } - - } - DUMMY_STATEMENT; - } break; - case 0x200C: /*UFS_IMAGE_DIMENSION_IN_BLOCK*/ {} break; - case 0x200D: /*UFS_IMAGE_BLOCK_HEADERS*/ {} break; - case 0x200E: /*UFS_IMAGE_BLOCK_COORDINATE*/ { - if (isyntax->parser.data_object_flags & (ISYNTAX_OBJECT_UFSImageBlockHeader | ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate)) { - isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; - // NOTE: the order of the coordinates stored in UFS_IMAGE_BLOCK_COORDINATE is variable and depends on UFS_IMAGE_DIMENSIONS_IN_CLUSTER (which might be stored later in the file) - // So, all we can do at this point is read them as 'raw coordinates' and link them with the actual coordinate once we have the complete UFS_IMAGE_DIMENSIONS_IN_CLUSTER info - parse_up_to_five_integers(value, (i32*)template->relative_coords_for_codeblock_in_cluster[isyntax->parser.block_header_index_for_cluster].raw_coords); - } - } break; - case 0x200F: /*UFS_IMAGE_BLOCK_COMPRESSION_METHOD*/ { - i32 compression_method = atoi(value); - if (compression_method == 16) { - image->compressor_version = 1; - } else if (compression_method == 19) { - image->compressor_version = 2; - } else { - // unknown compressor version - success = false; - } - } break; - case 0x2012: /*UFS_IMAGE_BLOCK_HEADER_TEMPLATE_ID*/ { - if (isyntax->parser.data_object_flags & (ISYNTAX_OBJECT_UFSImageBlockHeader | ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate)) { - isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; - isyntax_cluster_relative_coords_t* relative_coords_for_codeblock_in_cluster = template->relative_coords_for_codeblock_in_cluster + isyntax->parser.block_header_index_for_cluster; - relative_coords_for_codeblock_in_cluster->block_header_template_id = atoi(value); - } - } break; - case 0x2013: /*UFS_IMAGE_PIXEL_TRANSFORMATION_METHOD*/ {} break; - case 0x2014: { /*UFS_IMAGE_BLOCK_HEADER_TABLE*/ // data model <100 - // NOTE: mutually exclusive with UFS_IMAGE_BLOCK_HEADERS (either one or the other must be present) - size_t decoded_capacity = value_len; - size_t decoded_len = 0; - i32 last_char = value[value_len-1]; + free(decoded); + } + } break; + case 0x1013: /*DP_COLOR_MANAGEMENT*/ {} break; + case 0x1014: /*DP_IMAGE_POST_PROCESSING*/ {} break; + case 0x1015: /*DP_SHARPNESS_GAIN_RGB24*/ {} break; + case 0x1016: /*DP_CLAHE_CLIP_LIMIT_Y16*/ {} break; + case 0x1017: /*DP_CLAHE_NR_BINS_Y16*/ {} break; + case 0x1018: /*DP_CLAHE_CONTEXT_DIMENSION_Y16*/ {} break; + case 0x1019: /*DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR*/ {} break; + case 0x101A: /*DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL*/ {} break; + case 0x101B: /*DP_WAVELET_QUANTIZER*/ {} break; + case 0x101C: /*DP_WAVELET_DEADZONE*/ {} break; + case 0x1025: /*UFS_IMAGE_OPP_EXTREME_VERTEX*/ {} break; // data model >= 100 + case 0x2000: /*UFS_IMAGE_GENERAL_HEADERS*/ {} break; + case 0x2001: /*UFS_IMAGE_NUMBER_OF_BLOCKS*/ { + // NOTE: the actual number of stored codeblocks may be lower than this number. + image->number_of_blocks = atoi(value); + } break; + case 0x2002: /*UFS_IMAGE_DIMENSIONS_OVER_BLOCK*/ {} break; + case 0x2003: /*UFS_IMAGE_DIMENSIONS*/ {} break; + case 0x2004: /*UFS_IMAGE_DIMENSION_NAME*/ {} break; + case 0x2005: /*UFS_IMAGE_DIMENSION_TYPE*/ {} break; + case 0x2006: /*UFS_IMAGE_DIMENSION_UNIT*/ {} break; + case 0x2007: /*UFS_IMAGE_DIMENSION_SCALE_FACTOR*/ { + float mpp = atof(value); + if (isyntax->parser.dimension_index == 0 /*x*/) { + isyntax->mpp_x = mpp; + isyntax->is_mpp_known = true; + } else if (isyntax->parser.dimension_index == 1 /*y*/) { + isyntax->mpp_y = mpp; + isyntax->is_mpp_known = true; + } + } break; + case 0x2008: /*UFS_IMAGE_DIMENSION_DISCRETE_VALUES_STRING*/ {} break; + case 0x2009: /*UFS_IMAGE_BLOCK_HEADER_TEMPLATES*/ {} break; + case 0x200A: /*UFS_IMAGE_DIMENSION_RANGES*/ {} break; + case 0x200B: /*UFS_IMAGE_DIMENSION_RANGE*/ { + isyntax_image_dimension_range_t range = {0}; + parse_three_integers(value, &range.start, &range.step, &range.end); + i32 step_nonzero = (range.step != 0) ? range.step : 1; + range.numsteps = ((range.end + range.step) - range.start) / step_nonzero; + if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate) { + isyntax_block_header_template_t* template = isyntax->block_header_templates + isyntax->parser.block_header_template_index; + switch(isyntax->parser.dimension_index) { + default: break; + case 0: template->block_width = range.numsteps; break; + case 1: template->block_height = range.numsteps; break; + case 2: template->color_component = range.start; break; + case 3: template->scale = range.start; break; + case 4: template->waveletcoeff = (range.start == 0) ? 1 : 3; break; + } + DUMMY_STATEMENT; + } else if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageGeneralHeader) { + switch(isyntax->parser.dimension_index) { + default: break; + case 0: { + image->offset_x = range.start; + image->width = range.numsteps; + } break; + case 1: { + image->offset_y = range.start; + image->height = range.numsteps; + } break; + case 2: break; // always 3 color channels ("Y" "Co" "Cg"), no need to check + case 3: { + image->level_count = range.numsteps; + image->max_scale = range.numsteps - 1; + int32_t level_padding = (PER_LEVEL_PADDING << range.numsteps) - PER_LEVEL_PADDING; + image->width -= 2 * level_padding; + image->height -= 2 * level_padding; + } break; + case 4: break; // always 4 wavelet coefficients ("LL" "LH" "HL" "HH"), no need to check + } + DUMMY_STATEMENT; + } else if (isyntax->parser.data_object_flags & ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate) { + // The dimension range information within a cluster header template (specifically, the range.start) + // seems to provide the 'base values' for the dimension offsets for each of the codeblocks in the cluster template. + // To get the correct coordinate info for each codeblock, we need to add these 'base values' + // to the codeblock-specific coordinate values from the UFS_IMAGE_BLOCK_COORDINATE tags. + isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; + switch(isyntax->parser.dimension_index) { + default: break; + case 0: template->base_x = range.start; break; + case 1: template->base_y = range.start; break; + case 2: template->base_color_component = range.start; break; + case 3: template->base_scale = range.start; break; + case 4: template->base_waveletcoeff = range.start; break; + } + + } + DUMMY_STATEMENT; + } break; + case 0x200C: /*UFS_IMAGE_DIMENSION_IN_BLOCK*/ {} break; + case 0x200D: /*UFS_IMAGE_BLOCK_HEADERS*/ {} break; + case 0x200E: /*UFS_IMAGE_BLOCK_COORDINATE*/ { + if (isyntax->parser.data_object_flags & (ISYNTAX_OBJECT_UFSImageBlockHeader | ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate)) { + isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; + // NOTE: the order of the coordinates stored in UFS_IMAGE_BLOCK_COORDINATE is variable and depends on UFS_IMAGE_DIMENSIONS_IN_CLUSTER (which might be stored later in the file) + // So, all we can do at this point is read them as 'raw coordinates' and link them with the actual coordinate once we have the complete UFS_IMAGE_DIMENSIONS_IN_CLUSTER info + parse_up_to_five_integers(value, (i32*)template->relative_coords_for_codeblock_in_cluster[isyntax->parser.block_header_index_for_cluster].raw_coords); + } + } break; + case 0x200F: /*UFS_IMAGE_BLOCK_COMPRESSION_METHOD*/ { + i32 compression_method = atoi(value); + if (compression_method == 16) { + image->compressor_version = 1; + } else if (compression_method == 19) { + image->compressor_version = 2; + } else { + // unknown compressor version + success = false; + } + } break; + case 0x2012: /*UFS_IMAGE_BLOCK_HEADER_TEMPLATE_ID*/ { + if (isyntax->parser.data_object_flags & (ISYNTAX_OBJECT_UFSImageBlockHeader | ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate)) { + isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; + isyntax_cluster_relative_coords_t* relative_coords_for_codeblock_in_cluster = template->relative_coords_for_codeblock_in_cluster + isyntax->parser.block_header_index_for_cluster; + relative_coords_for_codeblock_in_cluster->block_header_template_id = atoi(value); + } + } break; + case 0x2013: /*UFS_IMAGE_PIXEL_TRANSFORMATION_METHOD*/ {} break; + case 0x2014: { /*UFS_IMAGE_BLOCK_HEADER_TABLE*/ // data model <100 + // NOTE: mutually exclusive with UFS_IMAGE_BLOCK_HEADERS (either one or the other must be present) + size_t decoded_capacity = value_len; + size_t decoded_len = 0; + i32 last_char = value[value_len-1]; /* FILE* test_out = file_stream_open_for_writing("test_b64.out"); file_stream_write(value, value_len, test_out); file_stream_close(test_out); */ - if (last_char == '/') { - value_len--; // The last character may cause the base64 decoding to fail if invalid - last_char = value[value_len-1]; - } - while (last_char == '\n' || last_char == '\r' || last_char == ' ') { - --value_len; - last_char = value[value_len-1]; - } - u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); - if (decoded) { - - u32 header_size = *(u32*) decoded + 0; - u8* block_header_start = decoded + 4; - isyntax_dicom_tag_header_t sequence_element = *(isyntax_dicom_tag_header_t*) (block_header_start); - if (sequence_element.size == 40) { - // We have a partial header structure, with 'Block Data Offset' and 'Block Size' missing (stored in Seektable) - // Full block header size (including the sequence element) is 48 bytes - u32 block_count = header_size / 48; - u32 should_be_zero = header_size % 48; - if (should_be_zero != 0) { - success = false; - } - - image->codeblock_count = block_count; - image->codeblocks = calloc(1, block_count * sizeof(isyntax_codeblock_t)); - image->header_codeblocks_are_partial = true; - - for (i32 i = 0; i < block_count; ++i) { - isyntax_partial_block_header_t* header = ((isyntax_partial_block_header_t*)(block_header_start)) + i; - isyntax_codeblock_t* codeblock = image->codeblocks + i; - codeblock->x_coordinate = header->x_coordinate; - codeblock->y_coordinate = header->y_coordinate; - codeblock->color_component = header->color_component; - codeblock->scale = header->scale; - codeblock->coefficient = header->coefficient; - codeblock->block_header_template_id = header->block_header_template_id; - DUMMY_STATEMENT; - } - - } else if (sequence_element.size == 72) { - // We have the complete header structure. (Nothing stored in Seektable) - u32 block_count = header_size / 80; - u32 should_be_zero = header_size % 80; - if (should_be_zero != 0) { - success = false; - } - - image->codeblock_count = block_count; - image->codeblocks = calloc(1, block_count * sizeof(isyntax_codeblock_t)); - image->header_codeblocks_are_partial = false; - - for (i32 i = 0; i < block_count; ++i) { - isyntax_full_block_header_t* header = ((isyntax_full_block_header_t*)(block_header_start)) + i; - isyntax_codeblock_t* codeblock = image->codeblocks + i; - codeblock->x_coordinate = header->x_coordinate; - codeblock->y_coordinate = header->y_coordinate; - codeblock->color_component = header->color_component; - codeblock->scale = header->scale; - codeblock->coefficient = header->coefficient; - codeblock->block_data_offset = header->block_data_offset; // extra - codeblock->block_size = header->block_size; // extra - codeblock->block_header_template_id = header->block_header_template_id; - DUMMY_STATEMENT; - } - } else { - success = false; - } - - free(decoded); - } else { - success = false; - } - } break; - case 0x2016: /*UFS_IMAGE_CLUSTER_HEADER_TEMPLATES*/ {} break; // data model >= 100 - case 0x2017: /*UFS_IMAGE_DIMENSIONS_OVER_CLUSTER*/ {} break; - case 0x201F: /*UFS_IMAGE_CLUSTER_HEADER_TABLE*/ { // data model >= 100 - size_t decoded_capacity = value_len; - size_t decoded_len = 0; - i32 last_char = value[value_len-1]; + if (last_char == '/') { + value_len--; // The last character may cause the base64 decoding to fail if invalid + last_char = value[value_len-1]; + } + while (last_char == '\n' || last_char == '\r' || last_char == ' ') { + --value_len; + last_char = value[value_len-1]; + } + u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); + if (decoded) { + + u32 header_size = *(u32*) decoded + 0; + u8* block_header_start = decoded + 4; + isyntax_dicom_tag_header_t sequence_element = *(isyntax_dicom_tag_header_t*) (block_header_start); + if (sequence_element.size == 40) { + // We have a partial header structure, with 'Block Data Offset' and 'Block Size' missing (stored in Seektable) + // Full block header size (including the sequence element) is 48 bytes + u32 block_count = header_size / 48; + u32 should_be_zero = header_size % 48; + if (should_be_zero != 0) { + success = false; + } + + image->codeblock_count = block_count; + image->codeblocks = calloc(1, block_count * sizeof(isyntax_codeblock_t)); + image->header_codeblocks_are_partial = true; + + for (i32 i = 0; i < block_count; ++i) { + isyntax_partial_block_header_t* header = ((isyntax_partial_block_header_t*)(block_header_start)) + i; + isyntax_codeblock_t* codeblock = image->codeblocks + i; + codeblock->x_coordinate = header->x_coordinate; + codeblock->y_coordinate = header->y_coordinate; + codeblock->color_component = header->color_component; + codeblock->scale = header->scale; + codeblock->coefficient = header->coefficient; + codeblock->block_header_template_id = header->block_header_template_id; + DUMMY_STATEMENT; + } + + } else if (sequence_element.size == 72) { + // We have the complete header structure. (Nothing stored in Seektable) + u32 block_count = header_size / 80; + u32 should_be_zero = header_size % 80; + if (should_be_zero != 0) { + success = false; + } + + image->codeblock_count = block_count; + image->codeblocks = calloc(1, block_count * sizeof(isyntax_codeblock_t)); + image->header_codeblocks_are_partial = false; + + for (i32 i = 0; i < block_count; ++i) { + isyntax_full_block_header_t* header = ((isyntax_full_block_header_t*)(block_header_start)) + i; + isyntax_codeblock_t* codeblock = image->codeblocks + i; + codeblock->x_coordinate = header->x_coordinate; + codeblock->y_coordinate = header->y_coordinate; + codeblock->color_component = header->color_component; + codeblock->scale = header->scale; + codeblock->coefficient = header->coefficient; + codeblock->block_data_offset = header->block_data_offset; // extra + codeblock->block_size = header->block_size; // extra + codeblock->block_header_template_id = header->block_header_template_id; + DUMMY_STATEMENT; + } + } else { + success = false; + } + + free(decoded); + } else { + success = false; + } + } break; + case 0x2016: /*UFS_IMAGE_CLUSTER_HEADER_TEMPLATES*/ {} break; // data model >= 100 + case 0x2017: /*UFS_IMAGE_DIMENSIONS_OVER_CLUSTER*/ {} break; + case 0x201F: /*UFS_IMAGE_CLUSTER_HEADER_TABLE*/ { // data model >= 100 + size_t decoded_capacity = value_len; + size_t decoded_len = 0; + i32 last_char = value[value_len-1]; /* FILE* test_out = file_stream_open_for_writing("test_b64.out"); file_stream_write(value, value_len, test_out); file_stream_close(test_out); */ - if (last_char == '/') { - value_len--; // The last character may cause the base64 decoding to fail if invalid - last_char = value[value_len-1]; - } - while (last_char == '\n' || last_char == '\r' || last_char == ' ') { - --value_len; - last_char = value[value_len-1]; - } - u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); - if (decoded) { - u8* decoded_end = decoded + decoded_len; - u32 header_size = *(u32*) decoded + 0; - u8* block_header_start = decoded + 4; - u8* pos = block_header_start; - isyntax_dicom_tag_header_t sequence_element = *(isyntax_dicom_tag_header_t*) pos; - - // first element should be a sequence tag - if (!(sequence_element.group == 0xfffe && sequence_element.element == 0xe000)) { - goto decoding_cluster_header_table_failed; - } - // pass 1: check how many clusters there are - i32 cluster_count = 0; - for (;;) { - u8* next_sequence_element_pos = pos + sizeof(isyntax_dicom_tag_header_t) + ((isyntax_dicom_tag_header_t*)pos)->size; - isyntax_dicom_tag_header_t* next_sequence_element = (isyntax_dicom_tag_header_t*)next_sequence_element_pos; - if (next_sequence_element_pos >= decoded_end || next_sequence_element->element != 0xe000) { - break; - } - ++cluster_count; - pos = next_sequence_element_pos; - } - - // preallocate memory for codeblocks and clusters - if (image->data_chunks == NULL) { - image->data_chunk_count = cluster_count; - image->data_chunks = calloc(1, cluster_count * sizeof(isyntax_data_chunk_t)); - } - ASSERT(image->number_of_blocks > 0); // TODO: handle error case - if (image->codeblocks == NULL) { - // NOTE: this value seems to be much larger than the actual number of codeblocks present in the file - image->codeblock_count = image->number_of_blocks; // from UFS_IMAGE_NUMBER_OF_BLOCKS attribute - image->codeblocks = calloc(1, image->codeblock_count * sizeof(isyntax_codeblock_t)); - } - - // pass 2: fill in all the information for each cluster - pos = block_header_start; - i32 running_codeblock_index = 0; - for (i32 i = 0; i < cluster_count; ++i) { - - sequence_element = *(isyntax_dicom_tag_header_t*) pos; - u8* next_sequence_element_pos = pos + sizeof(isyntax_dicom_tag_header_t) + ((isyntax_dicom_tag_header_t*)pos)->size; + if (last_char == '/') { + value_len--; // The last character may cause the base64 decoding to fail if invalid + last_char = value[value_len-1]; + } + while (last_char == '\n' || last_char == '\r' || last_char == ' ') { + --value_len; + last_char = value[value_len-1]; + } + u8* decoded = base64_decode((u8*)value, value_len, &decoded_len); + if (decoded) { + u8* decoded_end = decoded + decoded_len; + u32 header_size = *(u32*) decoded + 0; + u8* block_header_start = decoded + 4; + u8* pos = block_header_start; + isyntax_dicom_tag_header_t sequence_element = *(isyntax_dicom_tag_header_t*) pos; + + // first element should be a sequence tag + if (!(sequence_element.group == 0xfffe && sequence_element.element == 0xe000)) { + goto decoding_cluster_header_table_failed; + } + // pass 1: check how many clusters there are + i32 cluster_count = 0; + for (;;) { + u8* next_sequence_element_pos = pos + sizeof(isyntax_dicom_tag_header_t) + ((isyntax_dicom_tag_header_t*)pos)->size; + isyntax_dicom_tag_header_t* next_sequence_element = (isyntax_dicom_tag_header_t*)next_sequence_element_pos; + if (next_sequence_element_pos >= decoded_end || next_sequence_element->element != 0xe000) { + break; + } + ++cluster_count; + pos = next_sequence_element_pos; + } + + // preallocate memory for codeblocks and clusters + if (image->data_chunks == NULL) { + image->data_chunk_count = cluster_count; + image->data_chunks = calloc(1, cluster_count * sizeof(isyntax_data_chunk_t)); + } + ASSERT(image->number_of_blocks > 0); // TODO: handle error case + if (image->codeblocks == NULL) { + // NOTE: this value seems to be much larger than the actual number of codeblocks present in the file + image->codeblock_count = image->number_of_blocks; // from UFS_IMAGE_NUMBER_OF_BLOCKS attribute + image->codeblocks = calloc(1, image->codeblock_count * sizeof(isyntax_codeblock_t)); + } + + // pass 2: fill in all the information for each cluster + pos = block_header_start; + i32 running_codeblock_index = 0; + for (i32 i = 0; i < cluster_count; ++i) { + + sequence_element = *(isyntax_dicom_tag_header_t*) pos; + u8* next_sequence_element_pos = pos + sizeof(isyntax_dicom_tag_header_t) + ((isyntax_dicom_tag_header_t*)pos)->size; // dicom_tag_header_t* next_sequence_element = (dicom_tag_header_t*)next_sequence_element_pos; - u32 cluster_block_size = sequence_element.size; - u8* cluster_block_end = pos + sizeof(isyntax_dicom_tag_header_t) + cluster_block_size; - if (cluster_block_end > decoded_end) { - goto decoding_cluster_header_table_failed; // prevent out-of-bounds reading - } - - // advance to cluster coordinates - pos += sizeof(isyntax_dicom_tag_header_t); - isyntax_dicom_tag_header_t element = *(isyntax_dicom_tag_header_t*) pos; - u8* next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end) { - goto decoding_cluster_header_table_failed; - } - - i32 cluster_coordinate_count = element.size / 4; - i32* coordinates = (i32*)(pos + sizeof(isyntax_dicom_tag_header_t)); - if (cluster_coordinate_count < 2) { - // Expect only X and Y coordinates (provides the X/Y coordinates of the cluster) - // Scale, coefficient and color component are not needed, these can be derived from the cluster header templates - // (even so, scale seems to be included: in the example files, cluster_coordinate_count is 3) - goto decoding_cluster_header_table_failed; - } - i32 cluster_x = coordinates[0]; - i32 cluster_y = coordinates[1]; - - // read cluster header template ID - pos = next_element; - element = *(isyntax_dicom_tag_header_t*) pos; - next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end || element.size != 4) { - goto decoding_cluster_header_table_failed; - } - u32 cluster_header_template_id = *(u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); - if (cluster_header_template_id >= isyntax->cluster_header_template_count) { - goto decoding_cluster_header_table_failed; - } - isyntax_cluster_header_template_t* cluster_header_template = isyntax->cluster_header_templates + cluster_header_template_id; - if (cluster_coordinate_count >= 3) { - ASSERT(cluster_header_template->base_scale == coordinates[2]); - } - - // read cluster data offset - pos = next_element; - element = *(isyntax_dicom_tag_header_t*) pos; - next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end || element.size != 8) { - goto decoding_cluster_header_table_failed; - } - u64 cluster_data_offset = *(u64*)(pos + sizeof(isyntax_dicom_tag_header_t)); - - // read cluster size - pos = next_element; - element = *(isyntax_dicom_tag_header_t*) pos; - next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end || element.size != 8) { - goto decoding_cluster_header_table_failed; - } - u64 cluster_size = *(u64*)(pos + sizeof(isyntax_dicom_tag_header_t)); - - // read cluster block data offsets - pos = next_element; - element = *(isyntax_dicom_tag_header_t*) pos; - next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end) { - goto decoding_cluster_header_table_failed; - } - u32 block_count = element.size / 4; - u32* cluster_block_data_offsets = (u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); - - // read cluster block sizes - pos = next_element; - element = *(isyntax_dicom_tag_header_t*) pos; - next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; - if (next_element > cluster_block_end || element.size / 4 != block_count) { - goto decoding_cluster_header_table_failed; - } - u32* cluster_block_sizes = (u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); - - i32 top_codeblock_index = running_codeblock_index; - bool has_ll = false; - i32 highest_scale = 0; - ASSERT(running_codeblock_index + block_count <= image->codeblock_count); - for (i32 j = 0; j < block_count; ++j) { - isyntax_codeblock_t* codeblock = image->codeblocks + running_codeblock_index; - isyntax_cluster_relative_coords_t* relative_codeblock_in_cluster_info = cluster_header_template->relative_coords_for_codeblock_in_cluster + j; - codeblock->x_coordinate = cluster_x + relative_codeblock_in_cluster_info->x; - codeblock->y_coordinate = cluster_y + relative_codeblock_in_cluster_info->y; - codeblock->color_component = relative_codeblock_in_cluster_info->color_component; - codeblock->scale = relative_codeblock_in_cluster_info->scale; - if (codeblock->scale > highest_scale) highest_scale = codeblock->scale; - // account for different wavelet coefficient encoding in iSyntax v2 / data model >= 100 - codeblock->coefficient = (relative_codeblock_in_cluster_info->waveletcoeff == 3) ? 0 : 1; - if (codeblock->coefficient == 0) has_ll = true; - codeblock->block_data_offset = cluster_data_offset + cluster_block_data_offsets[j]; - codeblock->block_size = cluster_block_sizes[j]; - codeblock->block_header_template_id = relative_codeblock_in_cluster_info->block_header_template_id; - ++running_codeblock_index; - } - - isyntax_data_chunk_t* cluster = image->data_chunks + i; - cluster->offset = cluster_data_offset + cluster_block_data_offsets[0]; - cluster->size = cluster_size; - cluster->top_codeblock_index = top_codeblock_index; - cluster->codeblock_count_per_color = block_count / 3; - cluster->scale = highest_scale; - ASSERT(cluster->codeblock_count_per_color == isyntax_get_chunk_codeblocks_per_color_for_level(highest_scale, has_ll)); - - pos = next_sequence_element_pos; - } - - // TODO: prevent allocating too much memory in the first place? - if (running_codeblock_index < image->codeblock_count) { - // release excess allocated memory (there are fewer codeblocks present in the file than was anticipated) - image->codeblock_count = running_codeblock_index; - isyntax_codeblock_t* shrunk = realloc(image->codeblocks, image->codeblock_count * sizeof(isyntax_codeblock_t)); - ASSERT(shrunk); - image->codeblocks = shrunk; - } - - if (false) { decoding_cluster_header_table_failed: - success = false; - } - free(decoded); - - } else { - // base64 decoding failed - success = false; - } - } break; - case 0x2021: /*UFS_IMAGE_DIMENSIONS_IN_CLUSTER*/ { - // NOTE: Philips' documentation says this should have element tag 1021, but the example files have 2021! - isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; - template->dimension_count = parse_up_to_five_integers(value, template->dimension_order); - } break; // data model >= 100 - case 0x2023: /*UFS_IMAGE_VALID_DATA_ENVELOPES*/ {} break; // data model >= 100 - case 0x2024: /*UFS_IMAGE_OPP_EXTREME_VERTICES*/ {} break; // data model >= 100 - case 0x2025: /*UFS_IMAGE_OPP_EXTREME_VERTEX*/ {} break; // data model >= 100 - case 0x2026: /*UFS_IMAGE_VALID_ENVELOPE_DIMENSIONS*/ {} break; // data model >= 100 - case 0x2027: /*UFS_IMAGE_DIMENSION_ORIGIN*/ {} break; // data model >= 100 - case 0x2029: /*UFS_IMAGE_PIXEL_TRANSFORM_METHOD*/ {} break; // data model >= 100 - } - } break; - } - return success; + u32 cluster_block_size = sequence_element.size; + u8* cluster_block_end = pos + sizeof(isyntax_dicom_tag_header_t) + cluster_block_size; + if (cluster_block_end > decoded_end) { + goto decoding_cluster_header_table_failed; // prevent out-of-bounds reading + } + + // advance to cluster coordinates + pos += sizeof(isyntax_dicom_tag_header_t); + isyntax_dicom_tag_header_t element = *(isyntax_dicom_tag_header_t*) pos; + u8* next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end) { + goto decoding_cluster_header_table_failed; + } + + i32 cluster_coordinate_count = element.size / 4; + i32* coordinates = (i32*)(pos + sizeof(isyntax_dicom_tag_header_t)); + if (cluster_coordinate_count < 2) { + // Expect only X and Y coordinates (provides the X/Y coordinates of the cluster) + // Scale, coefficient and color component are not needed, these can be derived from the cluster header templates + // (even so, scale seems to be included: in the example files, cluster_coordinate_count is 3) + goto decoding_cluster_header_table_failed; + } + i32 cluster_x = coordinates[0]; + i32 cluster_y = coordinates[1]; + + // read cluster header template ID + pos = next_element; + element = *(isyntax_dicom_tag_header_t*) pos; + next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end || element.size != 4) { + goto decoding_cluster_header_table_failed; + } + u32 cluster_header_template_id = *(u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); + if (cluster_header_template_id >= isyntax->cluster_header_template_count) { + goto decoding_cluster_header_table_failed; + } + isyntax_cluster_header_template_t* cluster_header_template = isyntax->cluster_header_templates + cluster_header_template_id; + if (cluster_coordinate_count >= 3) { + ASSERT(cluster_header_template->base_scale == coordinates[2]); + } + + // read cluster data offset + pos = next_element; + element = *(isyntax_dicom_tag_header_t*) pos; + next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end || element.size != 8) { + goto decoding_cluster_header_table_failed; + } + u64 cluster_data_offset = *(u64*)(pos + sizeof(isyntax_dicom_tag_header_t)); + + // read cluster size + pos = next_element; + element = *(isyntax_dicom_tag_header_t*) pos; + next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end || element.size != 8) { + goto decoding_cluster_header_table_failed; + } + u64 cluster_size = *(u64*)(pos + sizeof(isyntax_dicom_tag_header_t)); + + // read cluster block data offsets + pos = next_element; + element = *(isyntax_dicom_tag_header_t*) pos; + next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end) { + goto decoding_cluster_header_table_failed; + } + u32 block_count = element.size / 4; + u32* cluster_block_data_offsets = (u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); + + // read cluster block sizes + pos = next_element; + element = *(isyntax_dicom_tag_header_t*) pos; + next_element = pos + sizeof(isyntax_dicom_tag_header_t) + element.size; + if (next_element > cluster_block_end || element.size / 4 != block_count) { + goto decoding_cluster_header_table_failed; + } + u32* cluster_block_sizes = (u32*)(pos + sizeof(isyntax_dicom_tag_header_t)); + + i32 top_codeblock_index = running_codeblock_index; + bool has_ll = false; + i32 highest_scale = 0; + ASSERT(running_codeblock_index + block_count <= image->codeblock_count); + for (i32 j = 0; j < block_count; ++j) { + isyntax_codeblock_t* codeblock = image->codeblocks + running_codeblock_index; + isyntax_cluster_relative_coords_t* relative_codeblock_in_cluster_info = cluster_header_template->relative_coords_for_codeblock_in_cluster + j; + codeblock->x_coordinate = cluster_x + relative_codeblock_in_cluster_info->x; + codeblock->y_coordinate = cluster_y + relative_codeblock_in_cluster_info->y; + codeblock->color_component = relative_codeblock_in_cluster_info->color_component; + codeblock->scale = relative_codeblock_in_cluster_info->scale; + if (codeblock->scale > highest_scale) highest_scale = codeblock->scale; + // account for different wavelet coefficient encoding in iSyntax v2 / data model >= 100 + codeblock->coefficient = (relative_codeblock_in_cluster_info->waveletcoeff == 3) ? 0 : 1; + if (codeblock->coefficient == 0) has_ll = true; + codeblock->block_data_offset = cluster_data_offset + cluster_block_data_offsets[j]; + codeblock->block_size = cluster_block_sizes[j]; + codeblock->block_header_template_id = relative_codeblock_in_cluster_info->block_header_template_id; + ++running_codeblock_index; + } + + isyntax_data_chunk_t* cluster = image->data_chunks + i; + cluster->offset = cluster_data_offset + cluster_block_data_offsets[0]; + cluster->size = cluster_size; + cluster->top_codeblock_index = top_codeblock_index; + cluster->codeblock_count_per_color = block_count / 3; + cluster->scale = highest_scale; + ASSERT(cluster->codeblock_count_per_color == isyntax_get_chunk_codeblocks_per_color_for_level(highest_scale, has_ll)); + + pos = next_sequence_element_pos; + } + + // TODO: prevent allocating too much memory in the first place? + if (running_codeblock_index < image->codeblock_count) { + // release excess allocated memory (there are fewer codeblocks present in the file than was anticipated) + image->codeblock_count = running_codeblock_index; + isyntax_codeblock_t* shrunk = realloc(image->codeblocks, image->codeblock_count * sizeof(isyntax_codeblock_t)); + ASSERT(shrunk); + image->codeblocks = shrunk; + } + + if (false) { decoding_cluster_header_table_failed: + success = false; + } + free(decoded); + + } else { + // base64 decoding failed + success = false; + } + } break; + case 0x2021: /*UFS_IMAGE_DIMENSIONS_IN_CLUSTER*/ { + // NOTE: Philips' documentation says this should have element tag 1021, but the example files have 2021! + isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + isyntax->parser.cluster_header_template_index; + template->dimension_count = parse_up_to_five_integers(value, template->dimension_order); + } break; // data model >= 100 + case 0x2023: /*UFS_IMAGE_VALID_DATA_ENVELOPES*/ {} break; // data model >= 100 + case 0x2024: /*UFS_IMAGE_OPP_EXTREME_VERTICES*/ {} break; // data model >= 100 + case 0x2025: /*UFS_IMAGE_OPP_EXTREME_VERTEX*/ {} break; // data model >= 100 + case 0x2026: /*UFS_IMAGE_VALID_ENVELOPE_DIMENSIONS*/ {} break; // data model >= 100 + case 0x2027: /*UFS_IMAGE_DIMENSION_ORIGIN*/ {} break; // data model >= 100 + case 0x2029: /*UFS_IMAGE_PIXEL_TRANSFORM_METHOD*/ {} break; // data model >= 100 + } + } break; + } + return success; } static bool validate_dicom_attr(const char* expected, const char* observed) { - bool ok = (strcmp(expected, observed) == 0); - if (!ok) { - console_print("iSyntax validation error: while reading DICOM metadata, expected '%s' but found '%s'\n", expected, observed); - } - return ok; + bool ok = (strcmp(expected, observed) == 0); + if (!ok) { + console_print("iSyntax validation error: while reading DICOM metadata, expected '%s' but found '%s'\n", expected, observed); + } + return ok; } void isyntax_xml_parser_init(isyntax_xml_parser_t* parser) { - parser->initialized = true; + parser->initialized = true; - parser->attrbuf_capacity = KILOBYTES(32); - parser->contentbuf_capacity = MEGABYTES(8); + parser->attrbuf_capacity = KILOBYTES(32); + parser->contentbuf_capacity = MEGABYTES(8); - parser->attrbuf = malloc(parser->attrbuf_capacity); - parser->attrbuf_end = parser->attrbuf + parser->attrbuf_capacity; - parser->attrcur = NULL; - parser->attrlen = 0; - parser->contentbuf = malloc(parser->contentbuf_capacity); - parser->contentcur = NULL; - parser->contentlen = 0; + parser->attrbuf = malloc(parser->attrbuf_capacity); + parser->attrbuf_end = parser->attrbuf + parser->attrbuf_capacity; + parser->attrcur = NULL; + parser->attrlen = 0; + parser->contentbuf = malloc(parser->contentbuf_capacity); + parser->contentcur = NULL; + parser->contentlen = 0; - parser->current_dicom_attribute_name[0] = '\0'; - parser->current_dicom_group_tag = 0; - parser->current_dicom_element_tag = 0; - parser->attribute_index = 0; - parser->current_node_type = ISYNTAX_NODE_NONE; + parser->current_dicom_attribute_name[0] = '\0'; + parser->current_dicom_group_tag = 0; + parser->current_dicom_element_tag = 0; + parser->attribute_index = 0; + parser->current_node_type = ISYNTAX_NODE_NONE; - // XML parsing using the yxml library. - // https://dev.yorhel.nl/yxml/man - size_t yxml_stack_buffer_size = KILOBYTES(32); - parser->x = (yxml_t*) malloc(sizeof(yxml_t) + yxml_stack_buffer_size); + // XML parsing using the yxml library. + // https://dev.yorhel.nl/yxml/man + size_t yxml_stack_buffer_size = KILOBYTES(32); + parser->x = (yxml_t*) malloc(sizeof(yxml_t) + yxml_stack_buffer_size); - yxml_init(parser->x, parser->x + 1, yxml_stack_buffer_size); + yxml_init(parser->x, parser->x + 1, yxml_stack_buffer_size); } static const char* get_spaces(i32 length) { - ASSERT(length >= 0); - static const char spaces[] = " "; - i32 spaces_len = COUNT(spaces) - 1; - i32 offset_from_end = MIN(spaces_len, length); - i32 offset = spaces_len - offset_from_end; - return spaces + offset; + ASSERT(length >= 0); + static const char spaces[] = " "; + i32 spaces_len = COUNT(spaces) - 1; + i32 offset_from_end = MIN(spaces_len, length); + i32 offset = spaces_len - offset_from_end; + return spaces + offset; } static void push_to_buffer_maybe_grow(u8** restrict dest, size_t* restrict dest_len, size_t* restrict dest_capacity, void* restrict src, size_t src_len) { - ASSERT(dest && dest_len && dest_capacity && src); - size_t old_len = *dest_len; - size_t new_len = old_len + src_len; - size_t capacity = *dest_capacity; - if (new_len > capacity) { - capacity = next_pow2(new_len); - u8* new_ptr = (u8*)realloc(*dest, capacity); - if (!new_ptr) panic(); - *dest = new_ptr; - *dest_capacity = capacity; - } - memcpy(*dest + old_len, src, src_len); - *dest_len = new_len; + ASSERT(dest && dest_len && dest_capacity && src); + size_t old_len = *dest_len; + size_t new_len = old_len + src_len; + size_t capacity = *dest_capacity; + if (new_len > capacity) { + capacity = next_pow2(new_len); + u8* new_ptr = (u8*)realloc(*dest, capacity); + if (!new_ptr) panic(); + *dest = new_ptr; + *dest_capacity = capacity; + } + memcpy(*dest + old_len, src, src_len); + *dest_len = new_len; } bool isyntax_parse_xml_header(isyntax_t* isyntax, char* xml_header, i64 chunk_length, bool is_last_chunk) { - yxml_t* x = NULL; - bool success = false; + yxml_t* x = NULL; + bool success = false; - static bool paranoid_mode = true; + static bool paranoid_mode = true; - isyntax_xml_parser_t* parser = &isyntax->parser; + isyntax_xml_parser_t* parser = &isyntax->parser; - if (!parser->initialized) { - isyntax_xml_parser_init(parser); - } - x = parser->x; + if (!parser->initialized) { + isyntax_xml_parser_init(parser); + } + x = parser->x; - if (0) { failed: cleanup: - if (parser->x) { - free(parser->x); - parser->x = NULL; - } - if (parser->attrbuf) { - free(parser->attrbuf); - parser->attrbuf = NULL; - } - if (parser->contentbuf) { - free(parser->contentbuf); - parser->contentbuf = NULL; - } - return success; - } + if (0) { failed: cleanup: + if (parser->x) { + free(parser->x); + parser->x = NULL; + } + if (parser->attrbuf) { + free(parser->attrbuf); + parser->attrbuf = NULL; + } + if (parser->contentbuf) { + free(parser->contentbuf); + parser->contentbuf = NULL; + } + return success; + } - // parse XML byte for byte - char* doc = xml_header; - for (i64 remaining_length = chunk_length; remaining_length > 0; --remaining_length, ++doc) { - int c = *doc; - if (c == '\0') { - // This should never trigger; iSyntax file is corrupt! - goto failed; - } - yxml_ret_t r = yxml_parse(x, c); - if (r == YXML_OK) { - continue; // nothing worthy of note has happened -> continue - } else if (r < 0) { - goto failed; - } else if (r > 0) { - // token - switch(r) { - case YXML_ELEMSTART: { - // start of an element: 'node_stack + parser->node_stack_index; - ++parser->node_stack_index; - isyntax_parser_node_t* node = parser->node_stack + parser->node_stack_index; - memset(node, 0, sizeof(isyntax_parser_node_t)); - // Inherit group and element of parent node - node->group = parent_node->group; - node->element = parent_node->element; - - parser->contentcur = parser->contentbuf; - *parser->contentcur = '\0'; - parser->contentlen = 0; - parser->attribute_index = 0; - if (strcmp(x->elem, "Attribute") == 0) { - node->node_type = ISYNTAX_NODE_LEAF; - } else if (strcmp(x->elem, "DataObject") == 0) { - node->node_type = ISYNTAX_NODE_BRANCH; - // push into the data object stack, to keep track of which type of DataObject we are parsing - // (we need this information to restore state when the XML element ends) - ++parser->data_object_stack_index; - parser->data_object_stack[parser->data_object_stack_index] = *parent_node; - // set relevant flag for which data object type we are now parsing - // NOTE: the data objects can have different DICOM groups, but currently there are no element - // ID collisions so we can switch on only the element ID. This may change in the future. - u32 flags = parser->data_object_flags; - switch(parent_node->element) { - default: break; - case 0: flags |= ISYNTAX_OBJECT_DPUfsImport; break; - case PIM_DP_SCANNED_IMAGES: flags |= ISYNTAX_OBJECT_DPScannedImage; break; - case UFS_IMAGE_GENERAL_HEADERS: flags |= ISYNTAX_OBJECT_UFSImageGeneralHeader; break; - case UFS_IMAGE_BLOCK_HEADER_TEMPLATES: flags |= ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate; break; - case UFS_IMAGE_DIMENSIONS: flags |= ISYNTAX_OBJECT_UFSImageDimension; break; - case UFS_IMAGE_DIMENSION_RANGES: flags |= ISYNTAX_OBJECT_UFSImageDimensionRange; break; - case DP_COLOR_MANAGEMENT: flags |= ISYNTAX_OBJECT_DPColorManagement; break; - case DP_IMAGE_POST_PROCESSING: flags |= ISYNTAX_OBJECT_DPImagePostProcessing; break; - case DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR: flags |= ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerColor; break; - case DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL: flags |= ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerLevel; break; - case PIIM_PIXEL_DATA_REPRESENTATION_SEQUENCE: flags |= ISYNTAX_OBJECT_PixelDataRepresentation; break; - case UFS_IMAGE_BLOCK_HEADERS: flags |= ISYNTAX_OBJECT_UFSImageBlockHeader; break; - case UFS_IMAGE_CLUSTER_HEADER_TEMPLATES: flags |= ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate; break; - case UFS_IMAGE_VALID_DATA_ENVELOPES: flags |= ISYNTAX_OBJECT_UFSImageValidDataEnvelope; break; - case UFS_IMAGE_OPP_EXTREME_VERTICES: flags |= ISYNTAX_OBJECT_UFSImageOppExtremeVertex; break; - } - parser->data_object_flags = flags; - } else if (strcmp(x->elem, "Array") == 0) { - node->node_type = ISYNTAX_NODE_ARRAY; - console_print_verbose("%sArray\n", get_spaces(parser->node_stack_index)); - } else { - node->node_type = ISYNTAX_NODE_NONE; - console_print_verbose("%selement start: %s\n", get_spaces(parser->node_stack_index), x->elem); - } - parser->current_node_type = node->node_type; - parser->current_node_has_children = false; - - } break; - - case YXML_CONTENT: { - // element content - if (!parser->contentcur) break; - - // Load iSyntax block header table (and other large XML tags) greedily and bypass yxml parsing overhead - if (parser->current_node_type == ISYNTAX_NODE_LEAF) { - u32 group = parser->current_dicom_group_tag; - u32 element = parser->current_dicom_element_tag; - isyntax_parser_node_t* node = parser->node_stack + parser->node_stack_index; - node->group = group; - node->element = element; - bool need_skip = (group == 0x301D && element == 0x2014) || // UFS_IMAGE_BLOCK_HEADER_TABLE - (group == 0x301D && element == 0x1005) || // PIM_DP_IMAGE_DATA - (group == 0x0028 && element == 0x2000); // DICOM_ICCPROFILE - - if (need_skip) { - parser->node_stack[parser->node_stack_index].has_base64_content = true; - char* content_start = doc; - char* pos = (char*)memchr(content_start, '<', remaining_length); - if (pos) { - i64 size = pos - content_start; - push_to_buffer_maybe_grow((u8**)&parser->contentbuf, &parser->contentlen, &parser->contentbuf_capacity, content_start, size); - parser->contentcur = parser->contentbuf + parser->contentlen; + // parse XML byte for byte + char* doc = xml_header; + for (i64 remaining_length = chunk_length; remaining_length > 0; --remaining_length, ++doc) { + int c = *doc; + if (c == '\0') { + // This should never trigger; iSyntax file is corrupt! + goto failed; + } + yxml_ret_t r = yxml_parse(x, c); + if (r == YXML_OK) { + continue; // nothing worthy of note has happened -> continue + } else if (r < 0) { + goto failed; + } else if (r > 0) { + // token + switch(r) { + case YXML_ELEMSTART: { + // start of an element: 'node_stack + parser->node_stack_index; + ++parser->node_stack_index; + isyntax_parser_node_t* node = parser->node_stack + parser->node_stack_index; + memset(node, 0, sizeof(isyntax_parser_node_t)); + // Inherit group and element of parent node + node->group = parent_node->group; + node->element = parent_node->element; + + parser->contentcur = parser->contentbuf; + *parser->contentcur = '\0'; + parser->contentlen = 0; + parser->attribute_index = 0; + if (strcmp(x->elem, "Attribute") == 0) { + node->node_type = ISYNTAX_NODE_LEAF; + } else if (strcmp(x->elem, "DataObject") == 0) { + node->node_type = ISYNTAX_NODE_BRANCH; + // push into the data object stack, to keep track of which type of DataObject we are parsing + // (we need this information to restore state when the XML element ends) + ++parser->data_object_stack_index; + parser->data_object_stack[parser->data_object_stack_index] = *parent_node; + // set relevant flag for which data object type we are now parsing + // NOTE: the data objects can have different DICOM groups, but currently there are no element + // ID collisions so we can switch on only the element ID. This may change in the future. + u32 flags = parser->data_object_flags; + switch(parent_node->element) { + default: break; + case 0: flags |= ISYNTAX_OBJECT_DPUfsImport; break; + case PIM_DP_SCANNED_IMAGES: flags |= ISYNTAX_OBJECT_DPScannedImage; break; + case UFS_IMAGE_GENERAL_HEADERS: flags |= ISYNTAX_OBJECT_UFSImageGeneralHeader; break; + case UFS_IMAGE_BLOCK_HEADER_TEMPLATES: flags |= ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate; break; + case UFS_IMAGE_DIMENSIONS: flags |= ISYNTAX_OBJECT_UFSImageDimension; break; + case UFS_IMAGE_DIMENSION_RANGES: flags |= ISYNTAX_OBJECT_UFSImageDimensionRange; break; + case DP_COLOR_MANAGEMENT: flags |= ISYNTAX_OBJECT_DPColorManagement; break; + case DP_IMAGE_POST_PROCESSING: flags |= ISYNTAX_OBJECT_DPImagePostProcessing; break; + case DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR: flags |= ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerColor; break; + case DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL: flags |= ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerLevel; break; + case PIIM_PIXEL_DATA_REPRESENTATION_SEQUENCE: flags |= ISYNTAX_OBJECT_PixelDataRepresentation; break; + case UFS_IMAGE_BLOCK_HEADERS: flags |= ISYNTAX_OBJECT_UFSImageBlockHeader; break; + case UFS_IMAGE_CLUSTER_HEADER_TEMPLATES: flags |= ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate; break; + case UFS_IMAGE_VALID_DATA_ENVELOPES: flags |= ISYNTAX_OBJECT_UFSImageValidDataEnvelope; break; + case UFS_IMAGE_OPP_EXTREME_VERTICES: flags |= ISYNTAX_OBJECT_UFSImageOppExtremeVertex; break; + } + parser->data_object_flags = flags; + } else if (strcmp(x->elem, "Array") == 0) { + node->node_type = ISYNTAX_NODE_ARRAY; + console_print_verbose("%sArray\n", get_spaces(parser->node_stack_index)); + } else { + node->node_type = ISYNTAX_NODE_NONE; + console_print_verbose("%selement start: %s\n", get_spaces(parser->node_stack_index), x->elem); + } + parser->current_node_type = node->node_type; + parser->current_node_has_children = false; + + } break; + + case YXML_CONTENT: { + // element content + if (!parser->contentcur) break; + + // Load iSyntax block header table (and other large XML tags) greedily and bypass yxml parsing overhead + if (parser->current_node_type == ISYNTAX_NODE_LEAF) { + u32 group = parser->current_dicom_group_tag; + u32 element = parser->current_dicom_element_tag; + isyntax_parser_node_t* node = parser->node_stack + parser->node_stack_index; + node->group = group; + node->element = element; + bool need_skip = (group == 0x301D && element == 0x2014) || // UFS_IMAGE_BLOCK_HEADER_TABLE + (group == 0x301D && element == 0x1005) || // PIM_DP_IMAGE_DATA + (group == 0x0028 && element == 0x2000); // DICOM_ICCPROFILE + + if (need_skip) { + parser->node_stack[parser->node_stack_index].has_base64_content = true; + char* content_start = doc; + char* pos = (char*)memchr(content_start, '<', remaining_length); + if (pos) { + i64 size = pos - content_start; + push_to_buffer_maybe_grow((u8**)&parser->contentbuf, &parser->contentlen, &parser->contentbuf_capacity, content_start, size); + parser->contentcur = parser->contentbuf + parser->contentlen; // console_print("iSyntax: skipped tag (0x%04x, 0x%04x) content length = %d\n", group, element, size); - doc += (size-1); // skip to the next tag - remaining_length -= (size-1); - break; - } else { + doc += (size-1); // skip to the next tag + remaining_length -= (size-1); + break; + } else { // console_print("iSyntax: skipped tag (0x%04x, 0x%04x) content length = %d\n", group, element, remaining_length); - push_to_buffer_maybe_grow((u8**)&parser->contentbuf, &parser->contentlen, &parser->contentbuf_capacity, content_start, remaining_length); - parser->contentcur = parser->contentbuf + parser->contentlen; - remaining_length = 0; // skip to the next chunk - break; - } - } - } - - char* tmp = x->data; - while (*tmp && parser->contentlen < parser->contentbuf_capacity) { - *(parser->contentcur++) = *(tmp++); - ++parser->contentlen; - // too long content -> resize buffer - if (parser->contentlen == parser->contentbuf_capacity) { - size_t new_capacity = parser->contentbuf_capacity * 2; - char* new_ptr = (char*)realloc(parser->contentbuf, new_capacity); - if (!new_ptr) panic(); - parser->contentbuf = new_ptr; - parser->contentcur = parser->contentbuf + parser->contentlen; - parser->contentbuf_capacity = new_capacity; + push_to_buffer_maybe_grow((u8**)&parser->contentbuf, &parser->contentlen, &parser->contentbuf_capacity, content_start, remaining_length); + parser->contentcur = parser->contentbuf + parser->contentlen; + remaining_length = 0; // skip to the next chunk + break; + } + } + } + + char* tmp = x->data; + while (*tmp && parser->contentlen < parser->contentbuf_capacity) { + *(parser->contentcur++) = *(tmp++); + ++parser->contentlen; + // too long content -> resize buffer + if (parser->contentlen == parser->contentbuf_capacity) { + size_t new_capacity = parser->contentbuf_capacity * 2; + char* new_ptr = (char*)realloc(parser->contentbuf, new_capacity); + if (!new_ptr) panic(); + parser->contentbuf = new_ptr; + parser->contentcur = parser->contentbuf + parser->contentlen; + parser->contentbuf_capacity = new_capacity; // console_print("isyntax_parse_xml_header(): XML content buffer overflow (resized buffer to %u)\n", new_capacity); - } - } + } + } - *parser->contentcur = '\0'; - } break; + *parser->contentcur = '\0'; + } break; - case YXML_ELEMEND: { - // end of an element: '.. />' or '' + case YXML_ELEMEND: { + // end of an element: '.. />' or '' - if (parser->current_node_type == ISYNTAX_NODE_LEAF && !parser->current_node_has_children) { - // Leaf node WITHOUT children. - // In this case we didn't already parse the attributes at the YXML_ATTREND stage. - // Now at the YXML_ELEMEND stage we can parse the complete tag at once (attributes + content). - console_print_verbose("%sDICOM: %-40s (0x%04x, 0x%04x), size:%-8u = %s\n", get_spaces(parser->node_stack_index), - parser->current_dicom_attribute_name, - parser->current_dicom_group_tag, parser->current_dicom_element_tag, - parser->contentlen, parser->contentbuf); + if (parser->current_node_type == ISYNTAX_NODE_LEAF && !parser->current_node_has_children) { + // Leaf node WITHOUT children. + // In this case we didn't already parse the attributes at the YXML_ATTREND stage. + // Now at the YXML_ELEMEND stage we can parse the complete tag at once (attributes + content). + console_print_verbose("%sDICOM: %-40s (0x%04x, 0x%04x), size:%-8u = %s\n", get_spaces(parser->node_stack_index), + parser->current_dicom_attribute_name, + parser->current_dicom_group_tag, parser->current_dicom_element_tag, + parser->contentlen, parser->contentbuf); // if (parser->node_stack[parser->node_stack_index].group == 0) { // DUMMY_STATEMENT; // probably the group is 0 because this is a top level node. // } - if (parser->node_stack_index == 2) { - isyntax_parse_ufsimport_child_node(isyntax, parser->current_dicom_group_tag, - parser->current_dicom_element_tag, - parser->contentbuf, parser->contentlen); - } else { - isyntax_parse_scannedimage_child_node(isyntax, parser->current_dicom_group_tag, - parser->current_dicom_element_tag, - parser->contentbuf, parser->contentlen); - } - } else { - // We have reached the end of a branch or array node, or a leaf node WITH children. - // In this case, the attributes have already been parsed at the YXML_ATTREND stage. - // Because their content is not text but child nodes, we do not need to touch the content buffer. - const char* elem_name = NULL; - if (parser->current_node_type == ISYNTAX_NODE_LEAF) { - // End of a leaf node WITH children. - elem_name = "Attribute"; - } else if (parser->current_node_type == ISYNTAX_NODE_BRANCH) { - elem_name = "DataObject"; - // pop data object stack - isyntax_parser_node_t data_object = parser->data_object_stack[parser->data_object_stack_index]; - --parser->data_object_stack_index; - // reset relevant data for the data object type we are now no longer parsing - u32 flags = parser->data_object_flags; - switch(data_object.element) { - default: break; - case 0: flags &= ~ISYNTAX_OBJECT_DPUfsImport; break; - case PIM_DP_SCANNED_IMAGES: { - flags &= ~ISYNTAX_OBJECT_DPScannedImage; - } break; - case UFS_IMAGE_GENERAL_HEADERS: { - flags &= ~ISYNTAX_OBJECT_UFSImageGeneralHeader; - parser->dimension_index = 0; - } break; - case UFS_IMAGE_BLOCK_HEADER_TEMPLATES: { - flags &= ~ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate; - ++parser->block_header_template_index; - ++isyntax->block_header_template_count; // TODO: refactor? - parser->dimension_index = 0; - } break; - case UFS_IMAGE_DIMENSIONS: { - flags &= ~ISYNTAX_OBJECT_UFSImageDimension; - ++parser->dimension_index; - } break; - case UFS_IMAGE_DIMENSION_RANGES: { - flags &= ~ISYNTAX_OBJECT_UFSImageDimensionRange; - ++parser->dimension_index; - } break; - case DP_COLOR_MANAGEMENT: flags &= ~ISYNTAX_OBJECT_DPColorManagement; break; - case DP_IMAGE_POST_PROCESSING: flags &= ~ISYNTAX_OBJECT_DPImagePostProcessing; break; - case DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR: flags &= ~ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerColor; break; - case DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL: flags &= ~ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerLevel; break; - case UFS_IMAGE_BLOCK_HEADERS: { - flags &= ~ISYNTAX_OBJECT_UFSImageBlockHeader; - if (flags & ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate) { - // NOTE: Within a UFSImageClusterHeaderTemplate, the UFSImageBlockHeader objects contain coordinate offset values that apply to that cluster template - // (Each UFSImageBlockHeader object corresponds to a codeblock in the cluster) - ++parser->block_header_index_for_cluster; - if (parser->block_header_index_for_cluster >= MAX_CODEBLOCKS_PER_CLUSTER) { - panic(); // TODO: unexpected error condition, fail more gracefully? - } - } - } break; - case UFS_IMAGE_CLUSTER_HEADER_TEMPLATES: { - // Finalize cluster header template object: fix up relative codeblock coordinates within cluster - // using the information from UFS_IMAGE_DIMENSION_RANGES (=base value) // TODO: verify this - // combined with the UFS_IMAGE_BLOCK_COORDINATE values (=offsets) using the dimension ordering info provided by UFS_IMAGE_DIMENSIONS_IN_CLUSTER - isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + parser->cluster_header_template_index; - template->codeblock_in_cluster_count = parser->block_header_index_for_cluster; - for (i32 i = 0; i < template->codeblock_in_cluster_count; ++i) { - isyntax_cluster_relative_coords_t* relative = template->relative_coords_for_codeblock_in_cluster + i; - // apply base values from UFS_IMAGE_DIMENSION_RANGES - relative->x = template->base_x; - relative->y = template->base_y; - relative->color_component = template->base_color_component; - relative->scale = template->base_scale; - relative->waveletcoeff = template->base_waveletcoeff; - // apply offset from UFS_IMAGE_BLOCK_COORDINATE - u32* dimensions_to_fix[5] = {&relative->x, &relative->y, &relative->color_component, &relative->scale, &relative->waveletcoeff}; - for (i32 dimension_index = 0; dimension_index < template->dimension_count; ++dimension_index) { - u32* to_fix = dimensions_to_fix[template->dimension_order[dimension_index]]; // get correct dimension from UFS_IMAGE_DIMENSIONS_IN_CLUSTER - *to_fix += relative->raw_coords[dimension_index]; // apply offset - } - } - // Pop flags and reset indices - flags &= ~ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate; - ++parser->cluster_header_template_index; - ++isyntax->cluster_header_template_count; // TODO: refactor? - parser->dimension_index = 0; - parser->block_header_index_for_cluster = 0; - } break; - case UFS_IMAGE_VALID_DATA_ENVELOPES: flags &= ~ISYNTAX_OBJECT_UFSImageValidDataEnvelope; break; - case UFS_IMAGE_OPP_EXTREME_VERTICES: flags &= ~ISYNTAX_OBJECT_UFSImageOppExtremeVertex; break; - } - parser->data_object_flags = flags; - } else if (parser->current_node_type == ISYNTAX_NODE_ARRAY) { - parser->dimension_index = 0; - elem_name = "Array"; - } - - console_print_verbose("%selement end: %s\n", get_spaces(parser->node_stack_index), elem_name); - } - - // 'Pop' context back to parent node - if (parser->node_stack_index > 0) { - --parser->node_stack_index; - parser->current_node_type = parser->node_stack[parser->node_stack_index].node_type; - parser->current_node_has_children = parser->node_stack[parser->node_stack_index].has_children; - } else { - //TODO: handle error condition - console_print_error("iSyntax XML error: closing element without matching start\n"); - } - - } break; - - case YXML_ATTRSTART: { - // attribute: 'Name=..' + if (parser->node_stack_index == 2) { + isyntax_parse_ufsimport_child_node(isyntax, parser->current_dicom_group_tag, + parser->current_dicom_element_tag, + parser->contentbuf, parser->contentlen); + } else { + isyntax_parse_scannedimage_child_node(isyntax, parser->current_dicom_group_tag, + parser->current_dicom_element_tag, + parser->contentbuf, parser->contentlen); + } + } else { + // We have reached the end of a branch or array node, or a leaf node WITH children. + // In this case, the attributes have already been parsed at the YXML_ATTREND stage. + // Because their content is not text but child nodes, we do not need to touch the content buffer. + const char* elem_name = NULL; + if (parser->current_node_type == ISYNTAX_NODE_LEAF) { + // End of a leaf node WITH children. + elem_name = "Attribute"; + } else if (parser->current_node_type == ISYNTAX_NODE_BRANCH) { + elem_name = "DataObject"; + // pop data object stack + isyntax_parser_node_t data_object = parser->data_object_stack[parser->data_object_stack_index]; + --parser->data_object_stack_index; + // reset relevant data for the data object type we are now no longer parsing + u32 flags = parser->data_object_flags; + switch(data_object.element) { + default: break; + case 0: flags &= ~ISYNTAX_OBJECT_DPUfsImport; break; + case PIM_DP_SCANNED_IMAGES: { + flags &= ~ISYNTAX_OBJECT_DPScannedImage; + } break; + case UFS_IMAGE_GENERAL_HEADERS: { + flags &= ~ISYNTAX_OBJECT_UFSImageGeneralHeader; + parser->dimension_index = 0; + } break; + case UFS_IMAGE_BLOCK_HEADER_TEMPLATES: { + flags &= ~ISYNTAX_OBJECT_UFSImageBlockHeaderTemplate; + ++parser->block_header_template_index; + ++isyntax->block_header_template_count; // TODO: refactor? + parser->dimension_index = 0; + } break; + case UFS_IMAGE_DIMENSIONS: { + flags &= ~ISYNTAX_OBJECT_UFSImageDimension; + ++parser->dimension_index; + } break; + case UFS_IMAGE_DIMENSION_RANGES: { + flags &= ~ISYNTAX_OBJECT_UFSImageDimensionRange; + ++parser->dimension_index; + } break; + case DP_COLOR_MANAGEMENT: flags &= ~ISYNTAX_OBJECT_DPColorManagement; break; + case DP_IMAGE_POST_PROCESSING: flags &= ~ISYNTAX_OBJECT_DPImagePostProcessing; break; + case DP_WAVELET_QUANTIZER_SETTINGS_PER_COLOR: flags &= ~ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerColor; break; + case DP_WAVELET_QUANTIZER_SETTINGS_PER_LEVEL: flags &= ~ISYNTAX_OBJECT_DPWaveletQuantizerSeetingsPerLevel; break; + case UFS_IMAGE_BLOCK_HEADERS: { + flags &= ~ISYNTAX_OBJECT_UFSImageBlockHeader; + if (flags & ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate) { + // NOTE: Within a UFSImageClusterHeaderTemplate, the UFSImageBlockHeader objects contain coordinate offset values that apply to that cluster template + // (Each UFSImageBlockHeader object corresponds to a codeblock in the cluster) + ++parser->block_header_index_for_cluster; + if (parser->block_header_index_for_cluster >= MAX_CODEBLOCKS_PER_CLUSTER) { + panic(); // TODO: unexpected error condition, fail more gracefully? + } + } + } break; + case UFS_IMAGE_CLUSTER_HEADER_TEMPLATES: { + // Finalize cluster header template object: fix up relative codeblock coordinates within cluster + // using the information from UFS_IMAGE_DIMENSION_RANGES (=base value) // TODO: verify this + // combined with the UFS_IMAGE_BLOCK_COORDINATE values (=offsets) using the dimension ordering info provided by UFS_IMAGE_DIMENSIONS_IN_CLUSTER + isyntax_cluster_header_template_t* template = isyntax->cluster_header_templates + parser->cluster_header_template_index; + template->codeblock_in_cluster_count = parser->block_header_index_for_cluster; + for (i32 i = 0; i < template->codeblock_in_cluster_count; ++i) { + isyntax_cluster_relative_coords_t* relative = template->relative_coords_for_codeblock_in_cluster + i; + // apply base values from UFS_IMAGE_DIMENSION_RANGES + relative->x = template->base_x; + relative->y = template->base_y; + relative->color_component = template->base_color_component; + relative->scale = template->base_scale; + relative->waveletcoeff = template->base_waveletcoeff; + // apply offset from UFS_IMAGE_BLOCK_COORDINATE + u32* dimensions_to_fix[5] = {&relative->x, &relative->y, &relative->color_component, &relative->scale, &relative->waveletcoeff}; + for (i32 dimension_index = 0; dimension_index < template->dimension_count; ++dimension_index) { + u32* to_fix = dimensions_to_fix[template->dimension_order[dimension_index]]; // get correct dimension from UFS_IMAGE_DIMENSIONS_IN_CLUSTER + *to_fix += relative->raw_coords[dimension_index]; // apply offset + } + } + // Pop flags and reset indices + flags &= ~ISYNTAX_OBJECT_UFSImageClusterHeaderTemplate; + ++parser->cluster_header_template_index; + ++isyntax->cluster_header_template_count; // TODO: refactor? + parser->dimension_index = 0; + parser->block_header_index_for_cluster = 0; + } break; + case UFS_IMAGE_VALID_DATA_ENVELOPES: flags &= ~ISYNTAX_OBJECT_UFSImageValidDataEnvelope; break; + case UFS_IMAGE_OPP_EXTREME_VERTICES: flags &= ~ISYNTAX_OBJECT_UFSImageOppExtremeVertex; break; + } + parser->data_object_flags = flags; + } else if (parser->current_node_type == ISYNTAX_NODE_ARRAY) { + parser->dimension_index = 0; + elem_name = "Array"; + } + + console_print_verbose("%selement end: %s\n", get_spaces(parser->node_stack_index), elem_name); + } + + // 'Pop' context back to parent node + if (parser->node_stack_index > 0) { + --parser->node_stack_index; + parser->current_node_type = parser->node_stack[parser->node_stack_index].node_type; + parser->current_node_has_children = parser->node_stack[parser->node_stack_index].has_children; + } else { + //TODO: handle error condition + console_print_error("iSyntax XML error: closing element without matching start\n"); + } + + } break; + + case YXML_ATTRSTART: { + // attribute: 'Name=..' // console_print_verbose("attr start: %s\n", x->attr); - parser->attrcur = parser->attrbuf; - *parser->attrcur = '\0'; - parser->attrlen = 0; - } break; - - case YXML_ATTRVAL: { - // attribute value - //console_print_verbose(" attr val: %s\n", x->attr); - if (!parser->attrcur) break; - char* tmp = x->data; - while (*tmp && parser->attrbuf < parser->attrbuf_end) { - *(parser->attrcur++) = *(tmp++); - ++parser->attrlen; - // too long content -> resize buffer - if (parser->attrlen == parser->attrbuf_capacity) { - size_t new_capacity = parser->attrbuf_capacity * 2; - char* new_ptr = (char*)realloc(parser->attrbuf, new_capacity); - if (!new_ptr) panic(); - parser->attrbuf = new_ptr; - parser->attrcur = parser->attrbuf + parser->attrlen; - parser->attrbuf_capacity = new_capacity; + parser->attrcur = parser->attrbuf; + *parser->attrcur = '\0'; + parser->attrlen = 0; + } break; + + case YXML_ATTRVAL: { + // attribute value + //console_print_verbose(" attr val: %s\n", x->attr); + if (!parser->attrcur) break; + char* tmp = x->data; + while (*tmp && parser->attrbuf < parser->attrbuf_end) { + *(parser->attrcur++) = *(tmp++); + ++parser->attrlen; + // too long content -> resize buffer + if (parser->attrlen == parser->attrbuf_capacity) { + size_t new_capacity = parser->attrbuf_capacity * 2; + char* new_ptr = (char*)realloc(parser->attrbuf, new_capacity); + if (!new_ptr) panic(); + parser->attrbuf = new_ptr; + parser->attrcur = parser->attrbuf + parser->attrlen; + parser->attrbuf_capacity = new_capacity; // console_print("isyntax_parse_xml_header(): XML attribute buffer overflow (resized buffer to %u)\n", new_capacity); - } - } - *parser->attrcur = '\0'; - } break; - - case YXML_ATTREND: { - // end of attribute '.."' - if (parser->attrcur) { - ASSERT(strlen(parser->attrbuf) == parser->attrlen); - - if (parser->current_node_type == ISYNTAX_NODE_LEAF) { - if (parser->attribute_index == 0 /* Name="..." */) { - if (paranoid_mode) validate_dicom_attr(x->attr, "Name"); - size_t copy_size = MIN(parser->attrlen, sizeof(parser->current_dicom_attribute_name)); - memcpy(parser->current_dicom_attribute_name, parser->attrbuf, copy_size); - i32 one_past_last_char = MIN(parser->attrlen, sizeof(parser->current_dicom_attribute_name)-1); - parser->current_dicom_attribute_name[one_past_last_char] = '\0'; - } else if (parser->attribute_index == 1 /* Group="0x...." */) { - if (paranoid_mode) validate_dicom_attr(x->attr, "Group"); - parser->current_dicom_group_tag = strtoul(parser->attrbuf, NULL, 0); - } else if (parser->attribute_index == 2 /* Element="0x...." */) { - if (paranoid_mode) validate_dicom_attr(x->attr, "Element"); - parser->current_dicom_element_tag = strtoul(parser->attrbuf, NULL, 0); - } else if (parser->attribute_index == 3 /* PMSVR="..." */) { - if (paranoid_mode) validate_dicom_attr(x->attr, "PMSVR"); - if (strcmp(parser->attrbuf, "IDataObjectArray") == 0) { - // Leaf node WITH children. - // Don't wait until YXML_ELEMEND to parse the attributes (this is the only opportunity we have!) - parser->current_node_has_children = true; - parser->node_stack[parser->node_stack_index].has_children = true; - console_print_verbose("%sDICOM: %-40s (0x%04x, 0x%04x), array\n", get_spaces(parser->node_stack_index), - parser->current_dicom_attribute_name, - parser->current_dicom_group_tag, parser->current_dicom_element_tag); - if (parser->node_stack_index == 2) { // At level of UfsImport - isyntax_parse_ufsimport_child_node(isyntax, parser->current_dicom_group_tag, - parser->current_dicom_element_tag, - parser->contentbuf, parser->contentlen); - } else { - bool parse_ok = isyntax_parse_scannedimage_child_node(isyntax, - parser->current_dicom_group_tag, - parser->current_dicom_element_tag, - parser->contentbuf, parser->contentlen); - } - - } - } - } else if (parser->current_node_type == ISYNTAX_NODE_BRANCH) { - // A DataObject node is supposed to have one attribute "ObjectType" - ASSERT(parser->attribute_index == 0); - ASSERT(strcmp(x->attr, "ObjectType") == 0); - console_print_verbose("%sDataObject %s = %s\n", get_spaces(parser->node_stack_index), x->attr, parser->attrbuf); - if (strcmp(parser->attrbuf, "DPScannedImage") == 0) { - // We started parsing a new image (which will be either a WSI, LABELIMAGE or MACROIMAGE). - parser->current_image = isyntax->images + isyntax->image_count; - parser->running_image_index = isyntax->image_count++; - } - } else { - console_print_verbose("%sattr %s = %s\n", get_spaces(parser->node_stack_index), x->attr, parser->attrbuf); - } - ++parser->attribute_index; - - } - } break; - - case YXML_PISTART: - case YXML_PICONTENT: - case YXML_PIEND: - break; // processing instructions (uninteresting, skip) - default: { - console_print_error("yxml_parse(): unrecognized token (%d)\n", r); - goto failed; - } - } - } - } + } + } + *parser->attrcur = '\0'; + } break; + + case YXML_ATTREND: { + // end of attribute '.."' + if (parser->attrcur) { + ASSERT(strlen(parser->attrbuf) == parser->attrlen); + + if (parser->current_node_type == ISYNTAX_NODE_LEAF) { + if (parser->attribute_index == 0 /* Name="..." */) { + if (paranoid_mode) validate_dicom_attr(x->attr, "Name"); + size_t copy_size = MIN(parser->attrlen, sizeof(parser->current_dicom_attribute_name)); + memcpy(parser->current_dicom_attribute_name, parser->attrbuf, copy_size); + i32 one_past_last_char = MIN(parser->attrlen, sizeof(parser->current_dicom_attribute_name)-1); + parser->current_dicom_attribute_name[one_past_last_char] = '\0'; + } else if (parser->attribute_index == 1 /* Group="0x...." */) { + if (paranoid_mode) validate_dicom_attr(x->attr, "Group"); + parser->current_dicom_group_tag = strtoul(parser->attrbuf, NULL, 0); + } else if (parser->attribute_index == 2 /* Element="0x...." */) { + if (paranoid_mode) validate_dicom_attr(x->attr, "Element"); + parser->current_dicom_element_tag = strtoul(parser->attrbuf, NULL, 0); + } else if (parser->attribute_index == 3 /* PMSVR="..." */) { + if (paranoid_mode) validate_dicom_attr(x->attr, "PMSVR"); + if (strcmp(parser->attrbuf, "IDataObjectArray") == 0) { + // Leaf node WITH children. + // Don't wait until YXML_ELEMEND to parse the attributes (this is the only opportunity we have!) + parser->current_node_has_children = true; + parser->node_stack[parser->node_stack_index].has_children = true; + console_print_verbose("%sDICOM: %-40s (0x%04x, 0x%04x), array\n", get_spaces(parser->node_stack_index), + parser->current_dicom_attribute_name, + parser->current_dicom_group_tag, parser->current_dicom_element_tag); + if (parser->node_stack_index == 2) { // At level of UfsImport + isyntax_parse_ufsimport_child_node(isyntax, parser->current_dicom_group_tag, + parser->current_dicom_element_tag, + parser->contentbuf, parser->contentlen); + } else { + bool parse_ok = isyntax_parse_scannedimage_child_node(isyntax, + parser->current_dicom_group_tag, + parser->current_dicom_element_tag, + parser->contentbuf, parser->contentlen); + } + + } + } + } else if (parser->current_node_type == ISYNTAX_NODE_BRANCH) { + // A DataObject node is supposed to have one attribute "ObjectType" + ASSERT(parser->attribute_index == 0); + ASSERT(strcmp(x->attr, "ObjectType") == 0); + console_print_verbose("%sDataObject %s = %s\n", get_spaces(parser->node_stack_index), x->attr, parser->attrbuf); + if (strcmp(parser->attrbuf, "DPScannedImage") == 0) { + // We started parsing a new image (which will be either a WSI, LABELIMAGE or MACROIMAGE). + parser->current_image = isyntax->images + isyntax->image_count; + parser->running_image_index = isyntax->image_count++; + } + } else { + console_print_verbose("%sattr %s = %s\n", get_spaces(parser->node_stack_index), x->attr, parser->attrbuf); + } + ++parser->attribute_index; + + } + } break; + + case YXML_PISTART: + case YXML_PICONTENT: + case YXML_PIEND: + break; // processing instructions (uninteresting, skip) + default: { + console_print_error("yxml_parse(): unrecognized token (%d)\n", r); + goto failed; + } + } + } + } - success = true; - if (is_last_chunk) { - goto cleanup; - } else { - return success; // no cleanup yet, we will still need resources until the last header chunk is reached. - } + success = true; + if (is_last_chunk) { + goto cleanup; + } else { + return success; // no cleanup yet, we will still need resources until the last header chunk is reached. + } } @@ -1217,22 +1220,22 @@ bool isyntax_parse_xml_header(isyntax_t* isyntax, char* xml_header, i64 chunk_le // https://stackoverflow.com/questions/21837008/how-to-convert-from-sign-magnitude-to-twos-complement // N.B. This function is its own inverse (conversion works the other way as well) static inline i16 signed_magnitude_to_twos_complement_16(u16 x) { - u16 m = -(x >> 15); - i16 result = (~m & x) | (((x & (u16)0x8000) - x) & m); - return result; + u16 m = -(x >> 15); + i16 result = (~m & x) | (((x & (u16)0x8000) - x) & m); + return result; } static inline i32 twos_complement_to_signed_magnitude(u32 x) { - u32 m = -(x >> 31); - i32 result = (~m & x) | (((x & 0x80000000) - x) & m); - return result; + u32 m = -(x >> 31); + i32 result = (~m & x) | (((x & 0x80000000) - x) & m); + return result; } static void signed_magnitude_to_twos_complement_16_block(u16* data, u32 len) { u32 aligned_len = (len / 8) * 8; u32 i = 0; #if defined(__SSE2__) - // Fast x86 SIMD version + // Fast x86 SIMD version for (; i < aligned_len; i += 8) { __m128i x = _mm_loadu_si128((__m128i*)(data + i)); __m128i sign_masks = _mm_srai_epi16(x, 15); // 0x0000 if positive, 0xFFFF if negative @@ -1255,10 +1258,10 @@ static void signed_magnitude_to_twos_complement_16_block(u16* data, u32 len) { } #endif // Slow version, for last unaligned elements or in case SIMD isn't available - for (; i < len; ++i) { - data[i] = signed_magnitude_to_twos_complement_16(data[i]); - } - ASSERT(i == len); + for (; i < len; ++i) { + data[i] = signed_magnitude_to_twos_complement_16(data[i]); + } + ASSERT(i == len); } // Convert a block of 16-bit signed integers to their absolute value @@ -1267,7 +1270,7 @@ static void signed_magnitude_to_absolute_value_16_block(i16* data, u32 len) { u32 aligned_len = (len / 8) * 8; u32 i = 0; #if defined(__SSE2__) - // Fast x86 SIMD version + // Fast x86 SIMD version for (; i < aligned_len; i += 8) { __m128i x = _mm_loadu_si128((__m128i*)(data + i)); __m128i sign_masks = _mm_srai_epi16(x, 15); // 0x0000 if positive, 0xFFFF if negative @@ -1292,10 +1295,10 @@ static void signed_magnitude_to_absolute_value_16_block(i16* data, u32 len) { } #endif // Slow version, for last unaligned elements or in case SIMD isn't available - for (; i < len; ++i) { - data[i] = (i16)(signed_magnitude_to_twos_complement_16((u16)data[i]) & 0x7FFF); - } - ASSERT(i == len); + for (; i < len; ++i) { + data[i] = (i16)(signed_magnitude_to_twos_complement_16((u16)data[i]) & 0x7FFF); + } + ASSERT(i == len); } @@ -1316,8 +1319,8 @@ void debug_convert_wavelet_coefficients_to_image2(icoeff_t* coefficients, i32 wi #if (DWT_COEFF_BITS==16) static u32 wavelet_coefficient_to_color_value(icoeff_t coefficient) { - u32 magnitude = ((u32)signed_magnitude_to_twos_complement_16(coefficient) & ~0x8000); - return magnitude; + u32 magnitude = ((u32)signed_magnitude_to_twos_complement_16(coefficient) & ~0x8000); + return magnitude; } #else static u32 wavelet_coefficient_to_color_value(icoeff_t coefficient) { @@ -1327,31 +1330,31 @@ static u32 wavelet_coefficient_to_color_value(icoeff_t coefficient) { #endif static rgba_t ycocg_to_rgb(i32 Y, i32 Co, i32 Cg) { - i32 tmp = Y - Cg/2; - i32 G = tmp + Cg; - i32 B = tmp - Co/2; - i32 R = B + Co; - return (rgba_t){ATMOST(255, R), ATMOST(255, G), ATMOST(255, B), 255}; + i32 tmp = Y - Cg/2; + i32 G = tmp + Cg; + i32 B = tmp - Co/2; + i32 R = B + Co; + return (rgba_t){ATMOST(255, R), ATMOST(255, G), ATMOST(255, B), 255}; } static rgba_t ycocg_to_bgr(i32 Y, i32 Co, i32 Cg) { - i32 tmp = Y - Cg/2; - i32 G = tmp + Cg; - i32 B = tmp - Co/2; - i32 R = B + Co; - return (rgba_t){ATMOST(255, B), ATMOST(255, G), ATMOST(255, R), 255}; + i32 tmp = Y - Cg/2; + i32 G = tmp + Cg; + i32 B = tmp - Co/2; + i32 R = B + Co; + return (rgba_t){ATMOST(255, B), ATMOST(255, G), ATMOST(255, R), 255}; } static u32* convert_ycocg_to_bgra_block(icoeff_t* Y, icoeff_t* Co, icoeff_t* Cg, i32 width, i32 height, i32 stride) { - i32 first_valid_pixel = ISYNTAX_IDWT_FIRST_VALID_PIXEL; - u32* bgra = (u32*)malloc(width * height * sizeof(u32)); // TODO: performance: block allocator + i32 first_valid_pixel = ISYNTAX_IDWT_FIRST_VALID_PIXEL; + u32* bgra = (u32*)malloc(width * height * sizeof(u32)); // TODO: performance: block allocator i32 aligned_width = (width / 8) * 8; - for (i32 y = 0; y < aligned_width; ++y) { - u32* dest = bgra + (y * width); + for (i32 y = 0; y < aligned_width; ++y) { + u32* dest = bgra + (y * width); i32 i = 0; #if defined(__SSE2__) && defined(__SSSE3__) - // Fast SIMD version (~2x faster on my system) + // Fast SIMD version (~2x faster on my system) for (; i < aligned_width; i += 8) { // Do the color space conversion __m128i Y_ = _mm_loadu_si128((__m128i*)(Y + i)); @@ -1404,15 +1407,15 @@ static u32* convert_ycocg_to_bgra_block(icoeff_t* Y, icoeff_t* Co, icoeff_t* Cg, } #endif // Slow version, for last unaligned elements or in case SIMD isn't available - for (; i < width; ++i) { - ((rgba_t*)dest)[i] = ycocg_to_bgr(Y[i], Co[i], Cg[i]); - } + for (; i < width; ++i) { + ((rgba_t*)dest)[i] = ycocg_to_bgr(Y[i], Co[i], Cg[i]); + } - Y += stride; - Co += stride; - Cg += stride; - } - return bgra; + Y += stride; + Co += stride; + Cg += stride; + } + return bgra; } static u32* convert_ycocg_to_rgba_block(icoeff_t* Y, icoeff_t* Co, icoeff_t* Cg, i32 width, i32 height, i32 stride) { @@ -1492,58 +1495,58 @@ static u32* convert_ycocg_to_rgba_block(icoeff_t* Y, icoeff_t* Co, icoeff_t* Cg, #define DEBUG_OUTPUT_IDWT_STEPS_AS_PNG 0 void isyntax_idwt(icoeff_t* idwt, i32 quadrant_width, i32 quadrant_height, bool output_steps_as_png, const char* png_name) { - i32 full_width = quadrant_width * 2; - i32 full_height= quadrant_height * 2; - i32 idwt_stride = full_width; + i32 full_width = quadrant_width * 2; + i32 full_height= quadrant_height * 2; + i32 idwt_stride = full_width; #if ISYNTAX_WANT_DEBUG_OUTPUT_PNG - if (output_steps_as_png) { + if (output_steps_as_png) { char filename[512]; snprintf(filename, sizeof(filename), "%s_step0.png", png_name); debug_convert_wavelet_coefficients_to_image2(idwt, full_width, full_height, filename); } #endif - // Horizontal pass - opj_dwt_t h = {0}; - size_t dwt_mem_size = (MAX(quadrant_width, quadrant_height)*2) * PARALLEL_COLS_53 * sizeof(icoeff_t); + // Horizontal pass + opj_dwt_t h = {0}; + size_t dwt_mem_size = (MAX(quadrant_width, quadrant_height)*2) * PARALLEL_COLS_53 * sizeof(icoeff_t); - h.mem = (icoeff_t*)alloca(dwt_mem_size); // TODO: need aligned memory? - h.sn = quadrant_width; // number of elements in low pass band - h.dn = quadrant_width; // number of elements in high pass band - h.cas = 1; + h.mem = (icoeff_t*)alloca(dwt_mem_size); // TODO: need aligned memory? + h.sn = quadrant_width; // number of elements in low pass band + h.dn = quadrant_width; // number of elements in high pass band + h.cas = 1; - for (i32 y = 0; y < full_height; ++y) { - icoeff_t* input_row = idwt + y * idwt_stride; - opj_idwt53_h(&h, input_row); - } + for (i32 y = 0; y < full_height; ++y) { + icoeff_t* input_row = idwt + y * idwt_stride; + opj_idwt53_h(&h, input_row); + } #if ISYNTAX_WANT_DEBUG_OUTPUT_PNG - if (output_steps_as_png) { + if (output_steps_as_png) { char filename[512]; snprintf(filename, sizeof(filename), "%s_step1.png", png_name); debug_convert_wavelet_coefficients_to_image2(idwt, full_width, full_height, filename); } #endif - // Vertical pass - opj_dwt_t v = {0}; - v.mem = h.mem; - v.sn = quadrant_height; // number of elements in low pass band - v.dn = quadrant_height; // number of elements in high pass band - v.cas = 1; - - i32 x; - i32 last_x = full_width; - for (x = 0; x + PARALLEL_COLS_53 <= last_x; x += PARALLEL_COLS_53) { - opj_idwt53_v(&v, idwt + x, idwt_stride, PARALLEL_COLS_53); - } - if (x < last_x) { - opj_idwt53_v(&v, idwt + x, idwt_stride, (last_x - x)); - } + // Vertical pass + opj_dwt_t v = {0}; + v.mem = h.mem; + v.sn = quadrant_height; // number of elements in low pass band + v.dn = quadrant_height; // number of elements in high pass band + v.cas = 1; + + i32 x; + i32 last_x = full_width; + for (x = 0; x + PARALLEL_COLS_53 <= last_x; x += PARALLEL_COLS_53) { + opj_idwt53_v(&v, idwt + x, idwt_stride, PARALLEL_COLS_53); + } + if (x < last_x) { + opj_idwt53_v(&v, idwt + x, idwt_stride, (last_x - x)); + } #if ISYNTAX_WANT_DEBUG_OUTPUT_PNG - if (output_steps_as_png) { + if (output_steps_as_png) { char filename[512]; snprintf(filename, sizeof(filename), "%s_step2.png", png_name); debug_convert_wavelet_coefficients_to_image2(idwt, full_width, full_height, filename); @@ -1553,497 +1556,497 @@ void isyntax_idwt(icoeff_t* idwt, i32 quadrant_width, i32 quadrant_height, bool } static inline void get_offsetted_coeff_blocks(icoeff_t** ll_hl_lh_hh, i32 offset, isyntax_tile_channel_t* color_channel, i32 block_stride, icoeff_t* black_dummy_coeff, icoeff_t* white_dummy_coeff) { - if (color_channel->coeff_ll) { - ll_hl_lh_hh[0] = color_channel->coeff_ll + offset; //ll - } else { - ll_hl_lh_hh[0] = white_dummy_coeff; - } - if (color_channel->coeff_h) { - ll_hl_lh_hh[1] = color_channel->coeff_h + offset; //hl - ll_hl_lh_hh[2] = color_channel->coeff_h + block_stride + offset; //lh - ll_hl_lh_hh[3] = color_channel->coeff_h + 2*block_stride + offset; //hh - } else { - ll_hl_lh_hh[1] = black_dummy_coeff; - ll_hl_lh_hh[2] = black_dummy_coeff; - ll_hl_lh_hh[3] = black_dummy_coeff; - } + if (color_channel->coeff_ll) { + ll_hl_lh_hh[0] = color_channel->coeff_ll + offset; //ll + } else { + ll_hl_lh_hh[0] = white_dummy_coeff; + } + if (color_channel->coeff_h) { + ll_hl_lh_hh[1] = color_channel->coeff_h + offset; //hl + ll_hl_lh_hh[2] = color_channel->coeff_h + block_stride + offset; //lh + ll_hl_lh_hh[3] = color_channel->coeff_h + 2*block_stride + offset; //hh + } else { + ll_hl_lh_hh[1] = black_dummy_coeff; + ll_hl_lh_hh[2] = black_dummy_coeff; + ll_hl_lh_hh[3] = black_dummy_coeff; + } } u32 isyntax_get_adjacent_tiles_mask(isyntax_level_t* level, i32 tile_x, i32 tile_y) { - ASSERT(tile_x >= 0 && tile_y >= 0); - ASSERT(tile_x < level->width_in_tiles && tile_y < level->height_in_tiles); - // 9 bits, corresponding to the surrounding tiles: - // 0x100 | 0x80 | 0x40 - // 0x20 | 0x10 | 8 - // 4 | 2 | 1 - u32 adj_tiles = 0x1FF; // all bits set - if (tile_y == 0) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_LEFT | ISYNTAX_ADJ_TILE_TOP_CENTER | ISYNTAX_ADJ_TILE_TOP_RIGHT); - if (tile_y == level->height_in_tiles-1) adj_tiles &= ~(ISYNTAX_ADJ_TILE_BOTTOM_LEFT | ISYNTAX_ADJ_TILE_BOTTOM_CENTER | ISYNTAX_ADJ_TILE_BOTTOM_RIGHT); - if (tile_x == 0) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_LEFT | ISYNTAX_ADJ_TILE_CENTER_LEFT | ISYNTAX_ADJ_TILE_BOTTOM_LEFT); - if (tile_x == level->width_in_tiles-1) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_RIGHT | ISYNTAX_ADJ_TILE_CENTER_RIGHT | ISYNTAX_ADJ_TILE_BOTTOM_RIGHT); - return adj_tiles; + ASSERT(tile_x >= 0 && tile_y >= 0); + ASSERT(tile_x < level->width_in_tiles && tile_y < level->height_in_tiles); + // 9 bits, corresponding to the surrounding tiles: + // 0x100 | 0x80 | 0x40 + // 0x20 | 0x10 | 8 + // 4 | 2 | 1 + u32 adj_tiles = 0x1FF; // all bits set + if (tile_y == 0) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_LEFT | ISYNTAX_ADJ_TILE_TOP_CENTER | ISYNTAX_ADJ_TILE_TOP_RIGHT); + if (tile_y == level->height_in_tiles-1) adj_tiles &= ~(ISYNTAX_ADJ_TILE_BOTTOM_LEFT | ISYNTAX_ADJ_TILE_BOTTOM_CENTER | ISYNTAX_ADJ_TILE_BOTTOM_RIGHT); + if (tile_x == 0) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_LEFT | ISYNTAX_ADJ_TILE_CENTER_LEFT | ISYNTAX_ADJ_TILE_BOTTOM_LEFT); + if (tile_x == level->width_in_tiles-1) adj_tiles &= ~(ISYNTAX_ADJ_TILE_TOP_RIGHT | ISYNTAX_ADJ_TILE_CENTER_RIGHT | ISYNTAX_ADJ_TILE_BOTTOM_RIGHT); + return adj_tiles; } u32 isyntax_get_adjacent_tiles_mask_only_existing(isyntax_level_t* level, i32 tile_x, i32 tile_y) { - u32 adjacent = isyntax_get_adjacent_tiles_mask(level, tile_x, tile_y); - u32 mask = 0; - if (adjacent & ISYNTAX_ADJ_TILE_TOP_LEFT) { - isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x-1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_LEFT; - } - if (adjacent & ISYNTAX_ADJ_TILE_TOP_CENTER) { - isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_CENTER; - } - if (adjacent & ISYNTAX_ADJ_TILE_TOP_RIGHT) { - isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x+1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_RIGHT; - } - if (adjacent & ISYNTAX_ADJ_TILE_CENTER_LEFT) { - isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x-1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER_LEFT; - } - if (adjacent & ISYNTAX_ADJ_TILE_CENTER) { - isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER; - } - if (adjacent & ISYNTAX_ADJ_TILE_CENTER_RIGHT) { - isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x+1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; - } - if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_LEFT) { - isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x-1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; - } - if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_CENTER) { - isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; - } - if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_RIGHT) { - isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x+1); - if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; - } - return mask; + u32 adjacent = isyntax_get_adjacent_tiles_mask(level, tile_x, tile_y); + u32 mask = 0; + if (adjacent & ISYNTAX_ADJ_TILE_TOP_LEFT) { + isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x-1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_LEFT; + } + if (adjacent & ISYNTAX_ADJ_TILE_TOP_CENTER) { + isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_CENTER; + } + if (adjacent & ISYNTAX_ADJ_TILE_TOP_RIGHT) { + isyntax_tile_t* tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x+1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_TOP_RIGHT; + } + if (adjacent & ISYNTAX_ADJ_TILE_CENTER_LEFT) { + isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x-1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER_LEFT; + } + if (adjacent & ISYNTAX_ADJ_TILE_CENTER) { + isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER; + } + if (adjacent & ISYNTAX_ADJ_TILE_CENTER_RIGHT) { + isyntax_tile_t* tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x+1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; + } + if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_LEFT) { + isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x-1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; + } + if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_CENTER) { + isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; + } + if (adjacent & ISYNTAX_ADJ_TILE_BOTTOM_RIGHT) { + isyntax_tile_t* tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x+1); + if (tile->exists) mask |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; + } + return mask; } u32 isyntax_idwt_tile_for_color_channel(isyntax_t* isyntax, isyntax_image_t* wsi, i32 scale, i32 tile_x, i32 tile_y, i32 color, icoeff_t* dest_buffer) { - isyntax_level_t* level = wsi->levels + scale; - ASSERT(tile_x >= 0 && tile_x < level->width_in_tiles); - ASSERT(tile_y >= 0 && tile_y < level->height_in_tiles); - isyntax_tile_t* tile = level->tiles + tile_y * level->width_in_tiles + tile_x; - isyntax_tile_channel_t* channel = tile->color_channels + color; + isyntax_level_t* level = wsi->levels + scale; + ASSERT(tile_x >= 0 && tile_x < level->width_in_tiles); + ASSERT(tile_y >= 0 && tile_y < level->height_in_tiles); + isyntax_tile_t* tile = level->tiles + tile_y * level->width_in_tiles + tile_x; + isyntax_tile_channel_t* channel = tile->color_channels + color; - u32 adj_tiles = isyntax_get_adjacent_tiles_mask(level, tile_x, tile_y); + u32 adj_tiles = isyntax_get_adjacent_tiles_mask(level, tile_x, tile_y); // ASSERT(channel->neighbors_loaded == adj_tiles); - // Prepare for stitching together the input image, with margins sampled from adjacent tiles for each quadrant - i32 pad_l = ISYNTAX_IDWT_PAD_L; - i32 pad_r = ISYNTAX_IDWT_PAD_R; - i32 pad_l_plus_r = pad_l + pad_r; + // Prepare for stitching together the input image, with margins sampled from adjacent tiles for each quadrant + i32 pad_l = ISYNTAX_IDWT_PAD_L; + i32 pad_r = ISYNTAX_IDWT_PAD_R; + i32 pad_l_plus_r = pad_l + pad_r; // ASSERT(sizeof(icoeff_t) * pad_amount == sizeof(u64)); // blit 64 bits == 4 pixels - i32 block_width = isyntax->block_width; - i32 block_height = isyntax->block_height; - i32 quadrant_width = block_width + pad_l_plus_r; - i32 quadrant_height = block_height + pad_l_plus_r; - i32 full_width = 2 * quadrant_width; - i32 full_height = 2 * quadrant_height; - icoeff_t* idwt = dest_buffer; // allocated/given by the caller ahead of time - - i32 dest_stride = full_width; - - // fill upper left quadrant with white - if (color == 0) { - for (i32 x = 0; x < quadrant_width; ++x) { - idwt[x] = 255; - } - for (i32 y = 1; y < quadrant_width; ++y) { - memcpy(idwt + y * dest_stride, idwt, quadrant_width * sizeof(icoeff_t)); - } - } - icoeff_t* h_dummy_coeff = isyntax->black_dummy_coeff; - icoeff_t* ll_dummy_coeff = (color == 0) ? isyntax->white_dummy_coeff : isyntax->black_dummy_coeff; - - - i32 source_stride = block_width; - i32 left_margin_source_x = block_width - pad_r; - i32 top_margin_source_y = block_height - pad_r; - size_t row_copy_size = block_width * sizeof(icoeff_t); - size_t pad_l_copy_size = pad_l * sizeof(icoeff_t); - size_t pad_r_copy_size = pad_r * sizeof(icoeff_t); - - i32 block_stride = block_width * block_height; - i32 quadrant_offsets[4] = {0, quadrant_width, full_width * quadrant_height, full_width * quadrant_height + quadrant_width}; - icoeff_t* quadrants[4] = {0}; - for (i32 i = 0; i < 4; ++i) { - quadrants[i] = idwt + quadrant_offsets[i]; - } + i32 block_width = isyntax->block_width; + i32 block_height = isyntax->block_height; + i32 quadrant_width = block_width + pad_l_plus_r; + i32 quadrant_height = block_height + pad_l_plus_r; + i32 full_width = 2 * quadrant_width; + i32 full_height = 2 * quadrant_height; + icoeff_t* idwt = dest_buffer; // allocated/given by the caller ahead of time + + i32 dest_stride = full_width; + + // fill upper left quadrant with white + if (color == 0) { + for (i32 x = 0; x < quadrant_width; ++x) { + idwt[x] = 255; + } + for (i32 y = 1; y < quadrant_width; ++y) { + memcpy(idwt + y * dest_stride, idwt, quadrant_width * sizeof(icoeff_t)); + } + } + icoeff_t* h_dummy_coeff = isyntax->black_dummy_coeff; + icoeff_t* ll_dummy_coeff = (color == 0) ? isyntax->white_dummy_coeff : isyntax->black_dummy_coeff; + + + i32 source_stride = block_width; + i32 left_margin_source_x = block_width - pad_r; + i32 top_margin_source_y = block_height - pad_r; + size_t row_copy_size = block_width * sizeof(icoeff_t); + size_t pad_l_copy_size = pad_l * sizeof(icoeff_t); + size_t pad_r_copy_size = pad_r * sizeof(icoeff_t); + + i32 block_stride = block_width * block_height; + i32 quadrant_offsets[4] = {0, quadrant_width, full_width * quadrant_height, full_width * quadrant_height + quadrant_width}; + icoeff_t* quadrants[4] = {0}; + for (i32 i = 0; i < 4; ++i) { + quadrants[i] = idwt + quadrant_offsets[i]; + } - icoeff_t* ll_hl_lh_hh[4] = {0}; - - u32 invalid_neighbors_ll = 0; - u32 invalid_neighbors_h = 0; - - // Now do the stitching, with margins sampled from adjacent tiles for each quadrant - // LL | HL - // LH | HH - - // top left corner - if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_LEFT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x-1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_LEFT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_LEFT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride) + left_margin_source_x, - color_channel, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t* source = ll_hl_lh_hh[i]; - icoeff_t* dest = quadrants[i]; - for (i32 y = 0; y < pad_l; ++y) { - memcpy(dest, source, pad_l_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } - // top center - if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_CENTER) { - isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + tile_x; - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_CENTER; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_CENTER; - get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride), - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + pad_l; - for (i32 y = 0; y < pad_l; ++y) { - memcpy(dest, source, row_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } - // top right corner - if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_RIGHT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x+1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_RIGHT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_RIGHT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride), - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + pad_l + block_width; - for (i32 y = 0; y < pad_l; ++y) { - memcpy(dest, source, pad_r_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } - // center left - if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER_LEFT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x-1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_CENTER_LEFT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_CENTER_LEFT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, left_margin_source_x, - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + (pad_l * dest_stride); - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, pad_l_copy_size); - dest += dest_stride; - source += source_stride; - } - } - } - } - // center (main tile) - if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER) { - get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, - channel, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + (pad_l * dest_stride) + pad_l; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, row_copy_size); - dest += dest_stride; - source += source_stride; - } - } + icoeff_t* ll_hl_lh_hh[4] = {0}; + + u32 invalid_neighbors_ll = 0; + u32 invalid_neighbors_h = 0; + + // Now do the stitching, with margins sampled from adjacent tiles for each quadrant + // LL | HL + // LH | HH + + // top left corner + if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_LEFT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x-1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_LEFT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_LEFT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride) + left_margin_source_x, + color_channel, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t* source = ll_hl_lh_hh[i]; + icoeff_t* dest = quadrants[i]; + for (i32 y = 0; y < pad_l; ++y) { + memcpy(dest, source, pad_l_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } + // top center + if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_CENTER) { + isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + tile_x; + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_CENTER; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_CENTER; + get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride), + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + pad_l; + for (i32 y = 0; y < pad_l; ++y) { + memcpy(dest, source, row_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } + // top right corner + if (adj_tiles & ISYNTAX_ADJ_TILE_TOP_RIGHT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y-1) * level->width_in_tiles + (tile_x+1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_TOP_RIGHT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_TOP_RIGHT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, (top_margin_source_y * source_stride), + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + pad_l + block_width; + for (i32 y = 0; y < pad_l; ++y) { + memcpy(dest, source, pad_r_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } + // center left + if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER_LEFT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x-1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_CENTER_LEFT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_CENTER_LEFT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, left_margin_source_x, + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + (pad_l * dest_stride); + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, pad_l_copy_size); + dest += dest_stride; + source += source_stride; + } + } + } + } + // center (main tile) + if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER) { + get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, + channel, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + (pad_l * dest_stride) + pad_l; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, row_copy_size); + dest += dest_stride; + source += source_stride; + } + } - } - // center right - if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER_RIGHT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x+1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + (pad_l * dest_stride) + pad_l + block_width; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, pad_r_copy_size); - dest += dest_stride; - source += source_stride; - } - } - } - } - // bottom left corner - if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_LEFT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x-1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, left_margin_source_x, - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride); - for (i32 y = 0; y < pad_r; ++y) { - memcpy(dest, source, pad_l_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } - // bottom center - if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_CENTER) { - isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + tile_x; - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; - get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride) + pad_l; - for (i32 y = 0; y < pad_r; ++y) { - memcpy(dest, source, row_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } - // bottom right corner - if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_RIGHT) { - isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x+1); - if (source_tile->exists) { - isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; - if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; - if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; - get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, - source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); - for (i32 i = 0; i < 4; ++i) { - icoeff_t *source = ll_hl_lh_hh[i]; - icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride) + pad_l + block_width; - for (i32 y = 0; y < pad_r; ++y) { - memcpy(dest, source, pad_r_copy_size); - source += source_stride; - dest += dest_stride; - } - } - } - } + } + // center right + if (adj_tiles & ISYNTAX_ADJ_TILE_CENTER_RIGHT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y) * level->width_in_tiles + (tile_x+1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_CENTER_RIGHT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + (pad_l * dest_stride) + pad_l + block_width; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, pad_r_copy_size); + dest += dest_stride; + source += source_stride; + } + } + } + } + // bottom left corner + if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_LEFT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x-1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_LEFT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, left_margin_source_x, + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride); + for (i32 y = 0; y < pad_r; ++y) { + memcpy(dest, source, pad_l_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } + // bottom center + if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_CENTER) { + isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + tile_x; + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_CENTER; + get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride) + pad_l; + for (i32 y = 0; y < pad_r; ++y) { + memcpy(dest, source, row_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } + // bottom right corner + if (adj_tiles & ISYNTAX_ADJ_TILE_BOTTOM_RIGHT) { + isyntax_tile_t* source_tile = level->tiles + (tile_y+1) * level->width_in_tiles + (tile_x+1); + if (source_tile->exists) { + isyntax_tile_channel_t* color_channel = source_tile->color_channels + color; + if (!color_channel->coeff_ll) invalid_neighbors_ll |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; + if (!color_channel->coeff_h) invalid_neighbors_h |= ISYNTAX_ADJ_TILE_BOTTOM_RIGHT; + get_offsetted_coeff_blocks(ll_hl_lh_hh, 0, + source_tile->color_channels + color, block_stride, h_dummy_coeff, ll_dummy_coeff); + for (i32 i = 0; i < 4; ++i) { + icoeff_t *source = ll_hl_lh_hh[i]; + icoeff_t *dest = quadrants[i] + ((pad_l + block_height) * dest_stride) + pad_l + block_width; + for (i32 y = 0; y < pad_r; ++y) { + memcpy(dest, source, pad_r_copy_size); + source += source_stride; + dest += dest_stride; + } + } + } + } - bool output_pngs = false; - const char* debug_png = "debug_idwt_"; - /* - if (scale == wsi->max_scale && tile_x == 1 && tile_y == 1 && color == 0) { - output_pngs = true; - }*/ - isyntax_idwt(idwt, quadrant_width, quadrant_height, output_pngs, debug_png); + bool output_pngs = false; + const char* debug_png = "debug_idwt_"; + /* + if (scale == wsi->max_scale && tile_x == 1 && tile_y == 1 && color == 0) { + output_pngs = true; + }*/ + isyntax_idwt(idwt, quadrant_width, quadrant_height, output_pngs, debug_png); - u32 invalid_edges = invalid_neighbors_h | invalid_neighbors_ll; - return invalid_edges; + u32 invalid_edges = invalid_neighbors_h | invalid_neighbors_ll; + return invalid_edges; } u32* isyntax_load_tile(isyntax_t* isyntax, isyntax_image_t* wsi, i32 scale, i32 tile_x, i32 tile_y, block_allocator_t* ll_coeff_block_allocator, bool decode_rgb) { - // printf("@@@ isyntax_load_tile scale=%d tile_x=%d tile_y=%d\n", scale, tile_x, tile_y); - isyntax_level_t* level = wsi->levels + scale; - ASSERT(tile_x >= 0 && tile_x < level->width_in_tiles); - ASSERT(tile_y >= 0 && tile_y < level->height_in_tiles); - isyntax_tile_t* tile = level->tiles + tile_y * level->width_in_tiles + tile_x; - i32 block_width = isyntax->block_width; - i32 block_height = isyntax->block_height; - size_t block_size = block_width * block_height * sizeof(icoeff_t); - i32 first_valid_pixel = ISYNTAX_IDWT_FIRST_VALID_PIXEL; - i32 idwt_width = 2 * (block_width + ISYNTAX_IDWT_PAD_L + ISYNTAX_IDWT_PAD_R); - i32 idwt_height = 2 * (block_height + ISYNTAX_IDWT_PAD_L + ISYNTAX_IDWT_PAD_R); - i32 idwt_stride = idwt_width; - size_t row_copy_size = block_width * sizeof(icoeff_t); - - temp_memory_t temp_memory = begin_temp_memory_on_local_thread(); - - icoeff_t* Y = NULL; - icoeff_t* Co = NULL; - icoeff_t* Cg = NULL; - - float elapsed_idwt = 0.0f; - float elapsed_malloc = 0.0f; - - u32 invalid_edges = 0; - - for (i32 color = 0; color < 3; ++color) { - i64 start_idwt = get_clock(); - // idwt will be allocated in temporary memory (only needed for the duration of this function) - size_t idwt_buffer_size = idwt_width * idwt_height * sizeof(icoeff_t); - icoeff_t* idwt = arena_push_size(temp_memory.arena, idwt_buffer_size); - memset(idwt, 0, idwt_buffer_size); - invalid_edges |= isyntax_idwt_tile_for_color_channel(isyntax, wsi, scale, tile_x, tile_y, color, idwt); - elapsed_idwt += get_seconds_elapsed(start_idwt, get_clock()); - ASSERT(idwt); - switch(color) { - case 0: Y = idwt; break; - case 1: Co = idwt; break; - case 2: Cg = idwt; break; - default: { - ASSERT(!"invalid code path"); - } break; - } + // printf("@@@ isyntax_load_tile scale=%d tile_x=%d tile_y=%d\n", scale, tile_x, tile_y); + isyntax_level_t* level = wsi->levels + scale; + ASSERT(tile_x >= 0 && tile_x < level->width_in_tiles); + ASSERT(tile_y >= 0 && tile_y < level->height_in_tiles); + isyntax_tile_t* tile = level->tiles + tile_y * level->width_in_tiles + tile_x; + i32 block_width = isyntax->block_width; + i32 block_height = isyntax->block_height; + size_t block_size = block_width * block_height * sizeof(icoeff_t); + i32 first_valid_pixel = ISYNTAX_IDWT_FIRST_VALID_PIXEL; + i32 idwt_width = 2 * (block_width + ISYNTAX_IDWT_PAD_L + ISYNTAX_IDWT_PAD_R); + i32 idwt_height = 2 * (block_height + ISYNTAX_IDWT_PAD_L + ISYNTAX_IDWT_PAD_R); + i32 idwt_stride = idwt_width; + size_t row_copy_size = block_width * sizeof(icoeff_t); + + temp_memory_t temp_memory = begin_temp_memory_on_local_thread(); + + icoeff_t* Y = NULL; + icoeff_t* Co = NULL; + icoeff_t* Cg = NULL; + + float elapsed_idwt = 0.0f; + float elapsed_malloc = 0.0f; + + u32 invalid_edges = 0; + + for (i32 color = 0; color < 3; ++color) { + i64 start_idwt = get_clock(); + // idwt will be allocated in temporary memory (only needed for the duration of this function) + size_t idwt_buffer_size = idwt_width * idwt_height * sizeof(icoeff_t); + icoeff_t* idwt = arena_push_size(temp_memory.arena, idwt_buffer_size); + memset(idwt, 0, idwt_buffer_size); + invalid_edges |= isyntax_idwt_tile_for_color_channel(isyntax, wsi, scale, tile_x, tile_y, color, idwt); + elapsed_idwt += get_seconds_elapsed(start_idwt, get_clock()); + ASSERT(idwt); + switch(color) { + case 0: Y = idwt; break; + case 1: Co = idwt; break; + case 2: Cg = idwt; break; + default: { + ASSERT(!"invalid code path"); + } break; + } - if (scale == 0) { - // No children to take care of at level 0. - continue; - } + if (scale == 0) { + // No children to take care of at level 0. + continue; + } - // Distribute result to child tiles if it was not distributed already. - isyntax_level_t* next_level = wsi->levels + (scale - 1); - isyntax_tile_t* child_top_left = next_level->tiles + (tile_y*2) * next_level->width_in_tiles + (tile_x*2); - isyntax_tile_t* child_top_right = child_top_left + 1; - isyntax_tile_t* child_bottom_left = child_top_left + next_level->width_in_tiles; - isyntax_tile_t* child_bottom_right = child_bottom_left + 1; + // Distribute result to child tiles if it was not distributed already. + isyntax_level_t* next_level = wsi->levels + (scale - 1); + isyntax_tile_t* child_top_left = next_level->tiles + (tile_y*2) * next_level->width_in_tiles + (tile_x*2); + isyntax_tile_t* child_top_right = child_top_left + 1; + isyntax_tile_t* child_bottom_left = child_top_left + next_level->width_in_tiles; + isyntax_tile_t* child_bottom_right = child_bottom_left + 1; - // TODO(avirodov): instead of releasing here, skip copy if still allocated. - if (child_top_left->color_channels[color].coeff_ll) { - block_free(ll_coeff_block_allocator, child_top_left->color_channels[color].coeff_ll); - } - if (child_top_right->color_channels[color].coeff_ll) { - block_free(ll_coeff_block_allocator, child_top_right->color_channels[color].coeff_ll); - } - if (child_bottom_left->color_channels[color].coeff_ll) { - block_free(ll_coeff_block_allocator, child_bottom_left->color_channels[color].coeff_ll); - } - if (child_bottom_right->color_channels[color].coeff_ll) { - block_free(ll_coeff_block_allocator, child_bottom_right->color_channels[color].coeff_ll); - } + // TODO(avirodov): instead of releasing here, skip copy if still allocated. + if (child_top_left->color_channels[color].coeff_ll) { + block_free(ll_coeff_block_allocator, child_top_left->color_channels[color].coeff_ll); + } + if (child_top_right->color_channels[color].coeff_ll) { + block_free(ll_coeff_block_allocator, child_top_right->color_channels[color].coeff_ll); + } + if (child_bottom_left->color_channels[color].coeff_ll) { + block_free(ll_coeff_block_allocator, child_bottom_left->color_channels[color].coeff_ll); + } + if (child_bottom_right->color_channels[color].coeff_ll) { + block_free(ll_coeff_block_allocator, child_bottom_right->color_channels[color].coeff_ll); + } - // NOTE: malloc() and free() can become a bottleneck, they don't scale well especially across many threads. - // We use a custom block allocator to address this. - i64 start_malloc = get_clock(); - child_top_left->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); - child_top_right->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); - child_bottom_left->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); - child_bottom_right->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); - elapsed_malloc += get_seconds_elapsed(start_malloc, get_clock()); - i32 dest_stride = block_width; - // Blit top left child LL block - { - icoeff_t* dest = child_top_left->color_channels[color].coeff_ll; - icoeff_t* source = idwt + (first_valid_pixel * idwt_stride) + first_valid_pixel; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, row_copy_size); - dest += dest_stride; - source += idwt_stride; - } - } - // Blit top right child LL block - { - icoeff_t* dest = child_top_right->color_channels[color].coeff_ll; - icoeff_t* source = idwt + (first_valid_pixel * idwt_stride) + first_valid_pixel + block_width; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, row_copy_size); - dest += dest_stride; - source += idwt_stride; - } - } - // Blit bottom left child LL block - { - icoeff_t* dest = child_bottom_left->color_channels[color].coeff_ll; - icoeff_t* source = idwt + ((first_valid_pixel + block_height) * idwt_stride) + first_valid_pixel; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, row_copy_size); - dest += dest_stride; - source += idwt_stride; - } - } - // Blit bottom right child LL block - { - icoeff_t* dest = child_bottom_right->color_channels[color].coeff_ll; - icoeff_t* source = idwt + ((first_valid_pixel + block_height) * idwt_stride) + first_valid_pixel + block_width; - for (i32 y = 0; y < block_height; ++y) { - memcpy(dest, source, row_copy_size); - dest += dest_stride; - source += idwt_stride; - } - } + // NOTE: malloc() and free() can become a bottleneck, they don't scale well especially across many threads. + // We use a custom block allocator to address this. + i64 start_malloc = get_clock(); + child_top_left->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); + child_top_right->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); + child_bottom_left->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); + child_bottom_right->color_channels[color].coeff_ll = (icoeff_t*)block_alloc(ll_coeff_block_allocator); + elapsed_malloc += get_seconds_elapsed(start_malloc, get_clock()); + i32 dest_stride = block_width; + // Blit top left child LL block + { + icoeff_t* dest = child_top_left->color_channels[color].coeff_ll; + icoeff_t* source = idwt + (first_valid_pixel * idwt_stride) + first_valid_pixel; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, row_copy_size); + dest += dest_stride; + source += idwt_stride; + } + } + // Blit top right child LL block + { + icoeff_t* dest = child_top_right->color_channels[color].coeff_ll; + icoeff_t* source = idwt + (first_valid_pixel * idwt_stride) + first_valid_pixel + block_width; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, row_copy_size); + dest += dest_stride; + source += idwt_stride; + } + } + // Blit bottom left child LL block + { + icoeff_t* dest = child_bottom_left->color_channels[color].coeff_ll; + icoeff_t* source = idwt + ((first_valid_pixel + block_height) * idwt_stride) + first_valid_pixel; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, row_copy_size); + dest += dest_stride; + source += idwt_stride; + } + } + // Blit bottom right child LL block + { + icoeff_t* dest = child_bottom_right->color_channels[color].coeff_ll; + icoeff_t* source = idwt + ((first_valid_pixel + block_height) * idwt_stride) + first_valid_pixel + block_width; + for (i32 y = 0; y < block_height; ++y) { + memcpy(dest, source, row_copy_size); + dest += dest_stride; + source += idwt_stride; + } + } - // After the last color channel, we can report that the children now have their LL blocks available. - if (color == 2) { - child_top_left->has_ll = true; - child_top_right->has_ll = true; - child_bottom_left->has_ll = true; - child_bottom_right->has_ll = true; - - if (invalid_edges != 0) { - console_print_error("load: scale=%d x=%d y=%d idwt time =%g invalid edges=%x\n", scale, tile_x, tile_y, elapsed_idwt, invalid_edges); - // early out - release_temp_memory(&temp_memory); - return NULL; - } - } - } + // After the last color channel, we can report that the children now have their LL blocks available. + if (color == 2) { + child_top_left->has_ll = true; + child_top_right->has_ll = true; + child_bottom_left->has_ll = true; + child_bottom_right->has_ll = true; + + if (invalid_edges != 0) { + console_print_error("load: scale=%d x=%d y=%d idwt time =%g invalid edges=%x\n", scale, tile_x, tile_y, elapsed_idwt, invalid_edges); + // early out + release_temp_memory(&temp_memory); + return NULL; + } + } + } - tile->is_loaded = true; // Meaning: it is now safe to start loading 'child' tiles of the next level - if (!decode_rgb) { - release_temp_memory(&temp_memory); // free Y, Co and Cg - return NULL; - } + tile->is_loaded = true; // Meaning: it is now safe to start loading 'child' tiles of the next level + if (!decode_rgb) { + release_temp_memory(&temp_memory); // free Y, Co and Cg + return NULL; + } - // For the Y (luminance) color channel, we actually need the absolute value of the Y-channel wavelet coefficient. - // (This doesn't hold for Co and Cg, those are are used directly as signed integers) + // For the Y (luminance) color channel, we actually need the absolute value of the Y-channel wavelet coefficient. + // (This doesn't hold for Co and Cg, those are are used directly as signed integers) signed_magnitude_to_absolute_value_16_block(Y, idwt_width * idwt_height); - // Reconstruct RGB image from separate color channels while cutting off margins - i64 start = get_clock(); - i32 tile_width = block_width * 2; - i32 tile_height = block_height * 2; + // Reconstruct RGB image from separate color channels while cutting off margins + i64 start = get_clock(); + i32 tile_width = block_width * 2; + i32 tile_height = block_height * 2; - i32 valid_offset = (first_valid_pixel * idwt_stride) + first_valid_pixel; - u32* bgra = convert_ycocg_to_bgra_block(Y + valid_offset, Co + valid_offset, Cg + valid_offset, - tile_width, tile_height, idwt_stride); - isyntax->total_rgb_transform_time += get_seconds_elapsed(start, get_clock()); + i32 valid_offset = (first_valid_pixel * idwt_stride) + first_valid_pixel; + u32* bgra = convert_ycocg_to_bgra_block(Y + valid_offset, Co + valid_offset, Cg + valid_offset, + tile_width, tile_height, idwt_stride); + isyntax->total_rgb_transform_time += get_seconds_elapsed(start, get_clock()); - // float elapsed_rgb = get_seconds_elapsed(start, get_clock()); - // console_print_verbose("load: scale=%d x=%d y=%d idwt time =%g rgb transform time=%g malloc time=%g\n", scale, tile_x, tile_y, elapsed_idwt, elapsed_rgb, elapsed_malloc); + // float elapsed_rgb = get_seconds_elapsed(start, get_clock()); + // console_print_verbose("load: scale=%d x=%d y=%d idwt time =%g rgb transform time=%g malloc time=%g\n", scale, tile_x, tile_y, elapsed_idwt, elapsed_rgb, elapsed_malloc); - /*if (scale == wsi->max_scale && tile_x == 1 && tile_y == 1) { - stbi_write_png("debug_dwt_output.png", tile_width, tile_height, 4, bgra, tile_width * 4); - }*/ + /*if (scale == wsi->max_scale && tile_x == 1 && tile_y == 1) { + stbi_write_png("debug_dwt_output.png", tile_width, tile_height, 4, bgra, tile_width * 4); + }*/ - release_temp_memory(&temp_memory); // free Y, Co and Cg - return bgra; + release_temp_memory(&temp_memory); // free Y, Co and Cg + return bgra; } @@ -2076,25 +2079,25 @@ u32* isyntax_load_tile(isyntax_t* isyntax, isyntax_image_t* wsi, i32 scale, i32 // The LL codeblock is only present at the highest scales. void isyntax_decompress_codeblock_in_chunk(isyntax_codeblock_t* codeblock, i32 block_width, i32 block_height, u8* chunk, u64 chunk_base_offset, i32 compressor_version, i16* out_buffer) { - i64 offset_in_chunk = codeblock->block_data_offset - chunk_base_offset; - ASSERT(offset_in_chunk >= 0); - isyntax_hulsken_decompress(chunk + offset_in_chunk, codeblock->block_size, - block_width, block_height, codeblock->coefficient, compressor_version, out_buffer); + i64 offset_in_chunk = codeblock->block_data_offset - chunk_base_offset; + ASSERT(offset_in_chunk >= 0); + isyntax_hulsken_decompress(chunk + offset_in_chunk, codeblock->block_size, + block_width, block_height, codeblock->coefficient, compressor_version, out_buffer); } // Read between 57 and 64 bits (7 bytes + 1-8 bits) from a bitstream (least significant bit first). // Requires that at least 7 safety bytes are present at the end of the stream (don't trigger a segmentation fault)! static inline u64 bitstream_lsb_read(u8* buffer, u32 pos) { - u64 raw = *(u64*)(buffer + pos / 8); - raw >>= pos % 8; - return raw; + u64 raw = *(u64*)(buffer + pos / 8); + raw >>= pos % 8; + return raw; } static inline u64 bitstream_lsb_read_advance(u8* buffer, i32* bits_read, i32 bits_to_read) { - u64 raw = *(u64*)(buffer + (*bits_read / 8)); - raw >>= (*bits_read / 8); - *bits_read += bits_to_read; - return raw; + u64 raw = *(u64*)(buffer + (*bits_read / 8)); + raw >>= (*bits_read / 8); + *bits_read += bits_to_read; + return raw; } // partly adapted from stb_image.h @@ -2104,22 +2107,22 @@ static inline u64 bitstream_lsb_read_advance(u8* buffer, i32* bits_read, i32 bit static const u16 size_bitmasks[17]={0,1,3,7,15,31,63,127,255,511,1023,2047,4095,8191,16383,32767,65535}; typedef struct huffman_t { - u16 fast[1 << HUFFMAN_FAST_BITS]; - u16 code[256]; - u8 size[256]; - u16 nonfast_symbols[256]; - u16 nonfast_code[256+7]; // extra safety bytes for SIMD vector operations - u16 nonfast_size[256]; - u16 nonfast_size_masks[256+7]; // extra safety bytes for SIMD vector operations + u16 fast[1 << HUFFMAN_FAST_BITS]; + u16 code[256]; + u8 size[256]; + u16 nonfast_symbols[256]; + u16 nonfast_code[256+7]; // extra safety bytes for SIMD vector operations + u16 nonfast_size[256]; + u16 nonfast_size_masks[256+7]; // extra safety bytes for SIMD vector operations } huffman_t; void save_code_in_huffman_fast_lookup_table(huffman_t* h, u32 code, u32 code_width, u8 symbol) { - ASSERT(code_width <= HUFFMAN_FAST_BITS); - i32 duplicate_bits = HUFFMAN_FAST_BITS - code_width; - for (u32 i = 0; i < (1 << duplicate_bits); ++i) { - u32 address = (i << code_width) | code; - h->fast[address] = symbol; - } + ASSERT(code_width <= HUFFMAN_FAST_BITS); + i32 duplicate_bits = HUFFMAN_FAST_BITS - code_width; + for (u32 i = 0; i < (1 << duplicate_bits); ++i) { + u32 address = (i << code_width) | code; + h->fast[address] = symbol; + } } //static u32 max_code_size; @@ -2164,228 +2167,228 @@ void isyntax_test_decompress_block(const char* filename) { */ bool isyntax_hulsken_decompress(u8* compressed, size_t compressed_size, i32 block_width, i32 block_height, - i32 coefficient, i32 compressor_version, i16* out_buffer) { - ASSERT(compressor_version == 1 || compressor_version == 2); - - // Read the header information stored in the codeblock. - // The layout varies depending on the version of the compressor used (version 1 or 2). - // All integers are stored little-endian, least-significant bit first. - // - // Version 1 layout: - // uint32 : serialized length (in bytes) - // uint8 : zero run symbol - // uint8 : zero run counter size (in bits) - // Version 2 layout: - // coeff_count (== 1 or 3) * coeff_bit_depth bits : 1 or 3 bitmasks, indicating which bitplanes are present - // uint8 : zero run symbol - // uint8 : zero run counter size (in bits) - // (variable length) : bitplane seektable (contains offsets to each of the bitplanes) - - // After the header section, the rest of the codeblock contains a Huffman tree, followed by a Huffman-coded - // message of 8-bit Huffman symbols, interspersed with 'zero run' symbols (for run-length encoding of zeroes). - - i32 coeff_count = (coefficient == 1) ? 3 : 1; - i32 coeff_bit_depth = 16; // fixed value for iSyntax - size_t coeff_buffer_size = coeff_count * block_width * block_height * sizeof(i16); - - // Early out if dummy/empty block - if (compressed_size <= 8) { - memset(out_buffer, 0, coeff_buffer_size); - return true; - } + i32 coefficient, i32 compressor_version, i16* out_buffer) { + ASSERT(compressor_version == 1 || compressor_version == 2); + + // Read the header information stored in the codeblock. + // The layout varies depending on the version of the compressor used (version 1 or 2). + // All integers are stored little-endian, least-significant bit first. + // + // Version 1 layout: + // uint32 : serialized length (in bytes) + // uint8 : zero run symbol + // uint8 : zero run counter size (in bits) + // Version 2 layout: + // coeff_count (== 1 or 3) * coeff_bit_depth bits : 1 or 3 bitmasks, indicating which bitplanes are present + // uint8 : zero run symbol + // uint8 : zero run counter size (in bits) + // (variable length) : bitplane seektable (contains offsets to each of the bitplanes) + + // After the header section, the rest of the codeblock contains a Huffman tree, followed by a Huffman-coded + // message of 8-bit Huffman symbols, interspersed with 'zero run' symbols (for run-length encoding of zeroes). + + i32 coeff_count = (coefficient == 1) ? 3 : 1; + i32 coeff_bit_depth = 16; // fixed value for iSyntax + size_t coeff_buffer_size = coeff_count * block_width * block_height * sizeof(i16); + + // Early out if dummy/empty block + if (compressed_size <= 8) { + memset(out_buffer, 0, coeff_buffer_size); + return true; + } - temp_memory_t temp_memory = begin_temp_memory_on_local_thread(); - - i32 bits_read = 0; - i32 block_size_in_bits = compressed_size * 8; - i64 serialized_length = 0; // In v1: stored in the first 4 bytes. In v2: derived calculation. - u32 bitmasks[3] = { 0x000FFFF, 0x000FFFF, 0x000FFFF }; // default in v1: all ones, can be overridden later - i32 total_mask_bits = coeff_bit_depth * coeff_count; - u8* byte_pos = compressed; - if (compressor_version == 1) { - serialized_length = *(u32*)byte_pos; - byte_pos += 4; - bits_read += 4*8; - } else { - if (coeff_count == 1) { - bitmasks[0] = *(u16*)(byte_pos); - byte_pos += 2; - bits_read += 2*8; - total_mask_bits = popcount(bitmasks[0]); - } else if (coeff_count == 3) { - bitmasks[0] = *(u16*)(byte_pos); - bitmasks[1] = *(u16*)(byte_pos+2); - bitmasks[2] = *(u16*)(byte_pos+4); - byte_pos += 6; - bits_read += 6*8; - total_mask_bits = popcount(bitmasks[0]) + popcount(bitmasks[1]) + popcount(bitmasks[2]); - } else { - panic("invalid coeff_count"); - } - serialized_length = total_mask_bits * (block_width * block_height / 8); - } + temp_memory_t temp_memory = begin_temp_memory_on_local_thread(); + + i32 bits_read = 0; + i32 block_size_in_bits = compressed_size * 8; + i64 serialized_length = 0; // In v1: stored in the first 4 bytes. In v2: derived calculation. + u32 bitmasks[3] = { 0x000FFFF, 0x000FFFF, 0x000FFFF }; // default in v1: all ones, can be overridden later + i32 total_mask_bits = coeff_bit_depth * coeff_count; + u8* byte_pos = compressed; + if (compressor_version == 1) { + serialized_length = *(u32*)byte_pos; + byte_pos += 4; + bits_read += 4*8; + } else { + if (coeff_count == 1) { + bitmasks[0] = *(u16*)(byte_pos); + byte_pos += 2; + bits_read += 2*8; + total_mask_bits = popcount(bitmasks[0]); + } else if (coeff_count == 3) { + bitmasks[0] = *(u16*)(byte_pos); + bitmasks[1] = *(u16*)(byte_pos+2); + bitmasks[2] = *(u16*)(byte_pos+4); + byte_pos += 6; + bits_read += 6*8; + total_mask_bits = popcount(bitmasks[0]) + popcount(bitmasks[1]) + popcount(bitmasks[2]); + } else { + panic("invalid coeff_count"); + } + serialized_length = total_mask_bits * (block_width * block_height / 8); + } - // Check that the serialized length is sane - if (serialized_length > 2 * coeff_buffer_size) { + // Check that the serialized length is sane + if (serialized_length > 2 * coeff_buffer_size) { // dump_block(compressed, compressed_size); - console_print_error("Error: isyntax_hulsken_decompress(): invalid codeblock, serialized_length too large (%lld)\n", serialized_length); - ASSERT(!"serialized_length too large"); - memset(out_buffer, 0, coeff_buffer_size); - release_temp_memory(&temp_memory); - return false; - } + console_print_error("Error: isyntax_hulsken_decompress(): invalid codeblock, serialized_length too large (%lld)\n", serialized_length); + ASSERT(!"serialized_length too large"); + memset(out_buffer, 0, coeff_buffer_size); + release_temp_memory(&temp_memory); + return false; + } - u8 zerorun_symbol = *(u8*)byte_pos++; - bits_read += 8; - u8 zero_counter_size = *(u8*)byte_pos++; - bits_read += 8; - - u32 bitplane_offsets[16] = {0}; - if (compressor_version >= 2) { - // Read bitplane seektable: a pointer is stored for a bit if it's represented in at least one of the bitmasks - u32 bitmasks_aggregate = 0; - for (i32 i = 0; i < coeff_count; ++i) { - bitmasks_aggregate |= bitmasks[i]; - } - i32 bitplane_ptr_count = popcount(bitmasks_aggregate); - i32 bitplane_ptr_bits = (i32)(log2f(serialized_length)) + 5; - u32 bitplane_ptr_mask = (1u << (bitplane_ptr_bits)) - 1; - for (i32 i = 0; i < bitplane_ptr_count - 1; ++i) { - u64 blob = bitstream_lsb_read(compressed, bits_read); - bitplane_offsets[i] = blob & bitplane_ptr_mask; - bits_read += bitplane_ptr_bits; - } - } + u8 zerorun_symbol = *(u8*)byte_pos++; + bits_read += 8; + u8 zero_counter_size = *(u8*)byte_pos++; + bits_read += 8; + + u32 bitplane_offsets[16] = {0}; + if (compressor_version >= 2) { + // Read bitplane seektable: a pointer is stored for a bit if it's represented in at least one of the bitmasks + u32 bitmasks_aggregate = 0; + for (i32 i = 0; i < coeff_count; ++i) { + bitmasks_aggregate |= bitmasks[i]; + } + i32 bitplane_ptr_count = popcount(bitmasks_aggregate); + i32 bitplane_ptr_bits = (i32)(log2f(serialized_length)) + 5; + u32 bitplane_ptr_mask = (1u << (bitplane_ptr_bits)) - 1; + for (i32 i = 0; i < bitplane_ptr_count - 1; ++i) { + u64 blob = bitstream_lsb_read(compressed, bits_read); + bitplane_offsets[i] = blob & bitplane_ptr_mask; + bits_read += bitplane_ptr_bits; + } + } - // Read Huffman table - huffman_t huffman = {0}; - memset(huffman.fast, 0x80, sizeof(huffman.fast)); - memset(huffman.nonfast_size_masks, 0xFF, sizeof(huffman.nonfast_size_masks)); - u32 fast_mask = (1 << HUFFMAN_FAST_BITS) - 1; - { - i32 code_size = 0; - u32 code = 0; - i32 nonfast_symbol_index = 0; - do { - if (bits_read >= block_size_in_bits) { + // Read Huffman table + huffman_t huffman = {0}; + memset(huffman.fast, 0x80, sizeof(huffman.fast)); + memset(huffman.nonfast_size_masks, 0xFF, sizeof(huffman.nonfast_size_masks)); + u32 fast_mask = (1 << HUFFMAN_FAST_BITS) - 1; + { + i32 code_size = 0; + u32 code = 0; + i32 nonfast_symbol_index = 0; + do { + if (bits_read >= block_size_in_bits) { // dump_block(compressed, compressed_size); - console_print_error("Error: isyntax_hulsken_decompress(): invalid codeblock, Huffman table extends out of bounds (compressed_size=%d)\n", compressed_size); - ASSERT(!"out of bounds"); - memset(out_buffer, 0, coeff_buffer_size); - release_temp_memory(&temp_memory); - return false; - } - // Read a chunk of bits large enough to 'always' have the whole Huffman code, followed by the 8-bit symbol. - // A blob of 57-64 bits is more than sufficient for a Huffman code of at most 16 bits. - // The bitstream is organized least significant bit first (treat as one giant little-endian integer). - // To read bits in the stream, look at the lowest bit positions. To advance the stream, shift right. - i32 bits_to_advance = 1; - u64 blob = bitstream_lsb_read(compressed, bits_read); // gives back between 57 and 64 bits. - - // 'Descend' into the tree until we hit a leaf node. - bool is_leaf = blob & 1; - // TODO: intrinsic? - while (!is_leaf) { - ++bits_to_advance; - blob >>= 1; - is_leaf = (blob & 1); - ++code_size; - } - blob >>= 1; - - // Read 8-bit Huffman symbol - u8 symbol = (u8)(blob); - huffman.code[symbol] = code; - huffman.size[symbol] = code_size; - - if (code_size <= HUFFMAN_FAST_BITS) { - // We can accelerate decoding of small Huffman codes by storing them in a lookup table. - // However, for the longer codes this becomes inefficient so in those cases we need another method. - save_code_in_huffman_fast_lookup_table(&huffman, code, code_size, symbol); + console_print_error("Error: isyntax_hulsken_decompress(): invalid codeblock, Huffman table extends out of bounds (compressed_size=%d)\n", compressed_size); + ASSERT(!"out of bounds"); + memset(out_buffer, 0, coeff_buffer_size); + release_temp_memory(&temp_memory); + return false; + } + // Read a chunk of bits large enough to 'always' have the whole Huffman code, followed by the 8-bit symbol. + // A blob of 57-64 bits is more than sufficient for a Huffman code of at most 16 bits. + // The bitstream is organized least significant bit first (treat as one giant little-endian integer). + // To read bits in the stream, look at the lowest bit positions. To advance the stream, shift right. + i32 bits_to_advance = 1; + u64 blob = bitstream_lsb_read(compressed, bits_read); // gives back between 57 and 64 bits. + + // 'Descend' into the tree until we hit a leaf node. + bool is_leaf = blob & 1; + // TODO: intrinsic? + while (!is_leaf) { + ++bits_to_advance; + blob >>= 1; + is_leaf = (blob & 1); + ++code_size; + } + blob >>= 1; + + // Read 8-bit Huffman symbol + u8 symbol = (u8)(blob); + huffman.code[symbol] = code; + huffman.size[symbol] = code_size; + + if (code_size <= HUFFMAN_FAST_BITS) { + // We can accelerate decoding of small Huffman codes by storing them in a lookup table. + // However, for the longer codes this becomes inefficient so in those cases we need another method. + save_code_in_huffman_fast_lookup_table(&huffman, code, code_size, symbol); // ++fast_count; - } else { - // Prepare the slow method for decoding Huffman codes that are too large to fit in the fast lookup table. - // (Slow method = iterating over possible symbols and checking if they match) - // Value being >= 256 in the 'fast' lookup table is the hint that we need the slow method. - u32 prefix = code & fast_mask; - u16 old_fast_data = huffman.fast[prefix]; - u8 old_lowest_symbol_index = old_fast_data & 0xFF; - u8 new_lowest_symbol_index = MIN(old_lowest_symbol_index, nonfast_symbol_index); - huffman.fast[prefix] = 256 + new_lowest_symbol_index; - huffman.nonfast_symbols[nonfast_symbol_index] = symbol; - huffman.nonfast_code[nonfast_symbol_index] = code; - huffman.nonfast_size[nonfast_symbol_index] = code_size; - huffman.nonfast_size_masks[nonfast_symbol_index] = size_bitmasks[code_size]; - ++nonfast_symbol_index; + } else { + // Prepare the slow method for decoding Huffman codes that are too large to fit in the fast lookup table. + // (Slow method = iterating over possible symbols and checking if they match) + // Value being >= 256 in the 'fast' lookup table is the hint that we need the slow method. + u32 prefix = code & fast_mask; + u16 old_fast_data = huffman.fast[prefix]; + u8 old_lowest_symbol_index = old_fast_data & 0xFF; + u8 new_lowest_symbol_index = MIN(old_lowest_symbol_index, nonfast_symbol_index); + huffman.fast[prefix] = 256 + new_lowest_symbol_index; + huffman.nonfast_symbols[nonfast_symbol_index] = symbol; + huffman.nonfast_code[nonfast_symbol_index] = code; + huffman.nonfast_size[nonfast_symbol_index] = code_size; + huffman.nonfast_size_masks[nonfast_symbol_index] = size_bitmasks[code_size]; + ++nonfast_symbol_index; // ++nonfast_count; - } - /*if (code_size > max_code_size) { - max_code_size = code_size; + } + /*if (code_size > max_code_size) { + max_code_size = code_size; // console_print("found the biggest code size: %d\n", code_size); - }*/ + }*/ // symbol_counts[symbol]++; - bits_to_advance += 8; - bits_read += bits_to_advance; - - // traverse back up the tree: find last zero -> flip to one - if (code_size == 0) { - break; // already done; this happens if there is only a root node, no leaves - } - u32 code_high_bit = (1 << (code_size - 1)); - bool found_zero = (~code) & code_high_bit; - while (!found_zero) { - --code_size; - if (code_size == 0) break; - code &= code_high_bit - 1; - code_high_bit >>= 1; - found_zero = (~code) & code_high_bit; - } - code |= code_high_bit; - } while(code_size > 0); - } + bits_to_advance += 8; + bits_read += bits_to_advance; + + // traverse back up the tree: find last zero -> flip to one + if (code_size == 0) { + break; // already done; this happens if there is only a root node, no leaves + } + u32 code_high_bit = (1 << (code_size - 1)); + bool found_zero = (~code) & code_high_bit; + while (!found_zero) { + --code_size; + if (code_size == 0) break; + code &= code_high_bit - 1; + code_high_bit >>= 1; + found_zero = (~code) & code_high_bit; + } + code |= code_high_bit; + } while(code_size > 0); + } - // Decode the message - u8* decompressed_buffer = (u8*)arena_push_size(temp_memory.arena, serialized_length); + // Decode the message + u8* decompressed_buffer = (u8*)arena_push_size(temp_memory.arena, serialized_length); - u32 zerorun_code = huffman.code[zerorun_symbol]; - u32 zerorun_code_size = huffman.size[zerorun_symbol]; - if (zerorun_code_size == 0) zerorun_code_size = 1; // handle special case of the 'empty' Huffman tree (root node is leaf node) - u32 zerorun_code_mask = (1 << zerorun_code_size) - 1; + u32 zerorun_code = huffman.code[zerorun_symbol]; + u32 zerorun_code_size = huffman.size[zerorun_symbol]; + if (zerorun_code_size == 0) zerorun_code_size = 1; // handle special case of the 'empty' Huffman tree (root node is leaf node) + u32 zerorun_code_mask = (1 << zerorun_code_size) - 1; - u32 zero_counter_mask = (1 << zero_counter_size) - 1; - i32 decompressed_length = 0; - while (bits_read < block_size_in_bits) { - if (decompressed_length >= serialized_length || bits_read >= block_size_in_bits) { - break; // done - } - i32 symbol = 0; - i32 code_size = 1; - u64 blob = bitstream_lsb_read(compressed, bits_read); - u32 fast_index = blob & fast_mask; - u16 c = huffman.fast[fast_index]; - if (c <= 255) { - // Lookup the symbol directly. - symbol = c; - code_size = huffman.size[symbol]; - } else { - bool match = false; - u8 lowest_possible_symbol_index = c & 0xFF; + u32 zero_counter_mask = (1 << zero_counter_size) - 1; + i32 decompressed_length = 0; + while (bits_read < block_size_in_bits) { + if (decompressed_length >= serialized_length || bits_read >= block_size_in_bits) { + break; // done + } + i32 symbol = 0; + i32 code_size = 1; + u64 blob = bitstream_lsb_read(compressed, bits_read); + u32 fast_index = blob & fast_mask; + u16 c = huffman.fast[fast_index]; + if (c <= 255) { + // Lookup the symbol directly. + symbol = c; + code_size = huffman.size[symbol]; + } else { + bool match = false; + u8 lowest_possible_symbol_index = c & 0xFF; #if !((defined(__SSE2__) && defined(__AVX__))) - for (i32 i = lowest_possible_symbol_index; i < 256; ++i) { - u8 test_size = huffman.nonfast_size[i]; - u16 test_code = huffman.nonfast_code[i]; - if ((blob & size_bitmasks[test_size]) == test_code) { - // match - code_size = test_size; - symbol = huffman.nonfast_symbols[i]; - match = true; - break; - } - } + for (i32 i = lowest_possible_symbol_index; i < 256; ++i) { + u8 test_size = huffman.nonfast_size[i]; + u16 test_code = huffman.nonfast_code[i]; + if ((blob & size_bitmasks[test_size]) == test_code) { + // match + code_size = test_size; + symbol = huffman.nonfast_symbols[i]; + match = true; + break; + } + } #else - // SIMD version using SSE2, very slightly faster than the version above + // SIMD version using SSE2, very slightly faster than the version above // (Need to compile with AVX enabled, otherwise unaligned loads will make it slower) // NOTE: Probably not a bottleneck, the loop below (nearly) always finishes after one iteration. for (i32 i = lowest_possible_symbol_index; i < 256; i += 8) { @@ -2407,204 +2410,204 @@ bool isyntax_hulsken_decompress(u8* compressed, size_t compressed_size, i32 bloc } DUMMY_STATEMENT; #endif - if (!match) { + if (!match) { // dump_block(compressed, compressed_size); - console_print_error("Error: isyntax_hulsken_decompress(): error decoding Huffman message (unknown symbol)\n"); - ASSERT(!"unknown symbol"); - memset(out_buffer, 0, coeff_buffer_size); - release_temp_memory(&temp_memory); - return false; - } - } - - if (code_size == 0) code_size = 1; // handle special case of the 'empty' Huffman tree (root node is leaf node) - - blob >>= code_size; - bits_read += code_size; - - // Handle run-length encoding of zeroes - if (symbol == zerorun_symbol) { - u32 numzeroes = blob & zero_counter_mask; - bits_read += zero_counter_size; - // A 'zero run' with length of zero means that this is not a zero run after all, but rather - // the 'escaped' zero run symbol itself which should be outputted. - if (numzeroes > 0) { - u32 actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one - if (decompressed_length + actual_numzeroes >= serialized_length || bits_read >= block_size_in_bits) { - // Reached the end, terminate - memset(decompressed_buffer + decompressed_length, 0, MIN(serialized_length - decompressed_length, actual_numzeroes)); - decompressed_length += actual_numzeroes; - break; - } - // If the next Huffman symbol is also the zero run symbol, then their counters actually refer to the same zero run. - // Basically, each extra zero run symbol expands the 'zero counter' bit depth, i.e.: - // n zero symbols -> depth becomes n * counter_bits - u32 total_zero_counter_size = zero_counter_size; - for(;;) { - // Peek ahead in the bitstream, grab any additional zero run symbols, and recalculate numzeroes. - blob = bitstream_lsb_read(compressed, bits_read); - u32 next_code = (blob & zerorun_code_mask); - if (next_code == zerorun_code) { - // The zero run continues - blob >>= zerorun_code_size; - u32 counter_extra_bits = blob & zero_counter_mask; - numzeroes <<= zero_counter_size; - numzeroes |= (counter_extra_bits); - total_zero_counter_size += zero_counter_size; - bits_read += zerorun_code_size + zero_counter_size; - actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one - if (decompressed_length + actual_numzeroes >= serialized_length || bits_read >= block_size_in_bits) { - break; // Reached the end, terminate - } - } else { - actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one - break; // no next zero run symbol, the zero run is finished - } - } + console_print_error("Error: isyntax_hulsken_decompress(): error decoding Huffman message (unknown symbol)\n"); + ASSERT(!"unknown symbol"); + memset(out_buffer, 0, coeff_buffer_size); + release_temp_memory(&temp_memory); + return false; + } + } - i32 bytes_to_write = MIN(serialized_length - decompressed_length, actual_numzeroes); - ASSERT(bytes_to_write > 0); - memset(decompressed_buffer + decompressed_length, 0, bytes_to_write); - decompressed_length += actual_numzeroes; - } else { - // This is not a 'zero run' after all, but an escaped symbol. So output the symbol. - decompressed_buffer[decompressed_length++] = symbol; - } - } else { - decompressed_buffer[decompressed_length++] = symbol; - } + if (code_size == 0) code_size = 1; // handle special case of the 'empty' Huffman tree (root node is leaf node) + + blob >>= code_size; + bits_read += code_size; + + // Handle run-length encoding of zeroes + if (symbol == zerorun_symbol) { + u32 numzeroes = blob & zero_counter_mask; + bits_read += zero_counter_size; + // A 'zero run' with length of zero means that this is not a zero run after all, but rather + // the 'escaped' zero run symbol itself which should be outputted. + if (numzeroes > 0) { + u32 actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one + if (decompressed_length + actual_numzeroes >= serialized_length || bits_read >= block_size_in_bits) { + // Reached the end, terminate + memset(decompressed_buffer + decompressed_length, 0, MIN(serialized_length - decompressed_length, actual_numzeroes)); + decompressed_length += actual_numzeroes; + break; + } + // If the next Huffman symbol is also the zero run symbol, then their counters actually refer to the same zero run. + // Basically, each extra zero run symbol expands the 'zero counter' bit depth, i.e.: + // n zero symbols -> depth becomes n * counter_bits + u32 total_zero_counter_size = zero_counter_size; + for(;;) { + // Peek ahead in the bitstream, grab any additional zero run symbols, and recalculate numzeroes. + blob = bitstream_lsb_read(compressed, bits_read); + u32 next_code = (blob & zerorun_code_mask); + if (next_code == zerorun_code) { + // The zero run continues + blob >>= zerorun_code_size; + u32 counter_extra_bits = blob & zero_counter_mask; + numzeroes <<= zero_counter_size; + numzeroes |= (counter_extra_bits); + total_zero_counter_size += zero_counter_size; + bits_read += zerorun_code_size + zero_counter_size; + actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one + if (decompressed_length + actual_numzeroes >= serialized_length || bits_read >= block_size_in_bits) { + break; // Reached the end, terminate + } + } else { + actual_numzeroes = (compressor_version == 2) ? numzeroes + 1 : numzeroes; // v2 stores actual count minus one + break; // no next zero run symbol, the zero run is finished + } + } + + i32 bytes_to_write = MIN(serialized_length - decompressed_length, actual_numzeroes); + ASSERT(bytes_to_write > 0); + memset(decompressed_buffer + decompressed_length, 0, bytes_to_write); + decompressed_length += actual_numzeroes; + } else { + // This is not a 'zero run' after all, but an escaped symbol. So output the symbol. + decompressed_buffer[decompressed_length++] = symbol; + } + } else { + decompressed_buffer[decompressed_length++] = symbol; + } - } + } - if (serialized_length != decompressed_length) { + if (serialized_length != decompressed_length) { // dump_block(compressed, compressed_size); - console_print("iSyntax: decompressed size mismatch (size=%d): expected %lld observed %d\n", - compressed_size, serialized_length, decompressed_length); - ASSERT(!"size mismatch"); - } - - i32 bytes_per_bitplane = (block_width * block_height) / 8; - if (compressor_version == 1) { - - i32 bytes_per_sample = 2; // ((coeff_bit_depth+7)/8); - i32 expected_bitmask_bits = (decompressed_length*8) / (block_width * block_height); + console_print("iSyntax: decompressed size mismatch (size=%d): expected %lld observed %d\n", + compressed_size, serialized_length, decompressed_length); + ASSERT(!"size mismatch"); + } - // Try to deduce the number of coefficients without knowing the header information - // NOTE: This is actually not necessary, because we do know coeff_count from the header. - i32 extra_bits = (decompressed_length*8) % (block_width * block_height); - if (extra_bits > 0) { - if (coeff_count != 1 && extra_bits == 1*16) { - coeff_count = 1; - } else if (coeff_count != 3 && extra_bits == 3*16) { - coeff_count = 3; - } - total_mask_bits = coeff_bit_depth * coeff_count; - } + i32 bytes_per_bitplane = (block_width * block_height) / 8; + if (compressor_version == 1) { + + i32 bytes_per_sample = 2; // ((coeff_bit_depth+7)/8); + i32 expected_bitmask_bits = (decompressed_length*8) / (block_width * block_height); + + // Try to deduce the number of coefficients without knowing the header information + // NOTE: This is actually not necessary, because we do know coeff_count from the header. + i32 extra_bits = (decompressed_length*8) % (block_width * block_height); + if (extra_bits > 0) { + if (coeff_count != 1 && extra_bits == 1*16) { + coeff_count = 1; + } else if (coeff_count != 3 && extra_bits == 3*16) { + coeff_count = 3; + } + total_mask_bits = coeff_bit_depth * coeff_count; + } - // If there are empty bitplanes: bitmasks stored at end of data - u64 expected_length = total_mask_bits * bytes_per_bitplane; - if (decompressed_length < expected_length) { - if (coeff_count == 1) { - bitmasks[0] = *(u16*)(decompressed_buffer + decompressed_length - 2); - total_mask_bits = popcount(bitmasks[0]); - } else if (coeff_count == 3) { - byte_pos = decompressed_buffer + decompressed_length - 6; - bitmasks[0] = *(u16*)(byte_pos); - bitmasks[1] = *(u16*)(byte_pos+2); - bitmasks[2] = *(u16*)(byte_pos+4); - total_mask_bits = popcount(bitmasks[0]) + popcount(bitmasks[1]) + popcount(bitmasks[2]); - } else { - panic("invalid coeff_count"); - } - expected_length = (total_mask_bits * block_width * block_height) / 8 + (coeff_count * 2); - ASSERT(decompressed_length == expected_length); - } - } + // If there are empty bitplanes: bitmasks stored at end of data + u64 expected_length = total_mask_bits * bytes_per_bitplane; + if (decompressed_length < expected_length) { + if (coeff_count == 1) { + bitmasks[0] = *(u16*)(decompressed_buffer + decompressed_length - 2); + total_mask_bits = popcount(bitmasks[0]); + } else if (coeff_count == 3) { + byte_pos = decompressed_buffer + decompressed_length - 6; + bitmasks[0] = *(u16*)(byte_pos); + bitmasks[1] = *(u16*)(byte_pos+2); + bitmasks[2] = *(u16*)(byte_pos+4); + total_mask_bits = popcount(bitmasks[0]) + popcount(bitmasks[1]) + popcount(bitmasks[2]); + } else { + panic("invalid coeff_count"); + } + expected_length = (total_mask_bits * block_width * block_height) / 8 + (coeff_count * 2); + ASSERT(decompressed_length == expected_length); + } + } - // unpack bitplanes - i32 compressed_bitplane_index = 0; - arena_align(temp_memory.arena, 32); - u16* coeff_buffer = (u16*)arena_push_size(temp_memory.arena, coeff_buffer_size); - memset(coeff_buffer, 0, coeff_buffer_size); - memset(out_buffer, 0, coeff_buffer_size); - - { - u32 running_bit_index = 0; - i32 running_coeff_index = 0; - u32 bitmasks_copy[3]; - memcpy(bitmasks_copy, bitmasks, sizeof(bitmasks)); - for (i32 bitplane_index = 0; bitplane_index < total_mask_bits; ++bitplane_index) { - u8* bitplane = decompressed_buffer + (bitplane_index * bytes_per_bitplane); - - // horribly complicated 'for loop'-style iteration, needed because v1 and v2 store bitplanes in a different order - // v1: iterate over the bitplanes (each coefficient separately) - // v2: for each bitplane, alternate the coefficients (striped) - // TODO: refactor? - if (compressor_version == 1) { - // find next bit for current coeff - for (;;) { - if (running_coeff_index >= coeff_count) { - panic("too many bitplanes"); - } - u16 bitmask = bitmasks_copy[running_coeff_index]; - if (bitmask) { - running_bit_index = bit_scan_forward(bitmask); - ASSERT(running_bit_index < 16); - bitmasks_copy[running_coeff_index] &= ~(1 << running_bit_index); - break; // success - } else { - ++running_coeff_index; - } - } - } else { - // compressor version 2: alternating coefficients - for (;;) { - if (running_bit_index >= 16) { - panic("too many bitplanes"); - } - if (running_coeff_index < coeff_count) { - u16 bitmask = bitmasks_copy[running_coeff_index]; - if (bitmask & (1 << running_bit_index)) { - bitmasks_copy[running_coeff_index] &= ~(1 << running_bit_index); - break; // success - } else { - ++running_coeff_index; // keep looking - } - } else { - running_coeff_index = 0; - ++running_bit_index; // keep looking - } - } - } - // Now we figured out which coeff and bit number this bitplane belongs to - - u16* current_coeff_buffer = coeff_buffer + (running_coeff_index * (block_width * block_height)); - u16* current_out_buffer = (u16*)out_buffer + (running_coeff_index * (block_width * block_height)); - - // Do the bitplane unpacking - for (i32 i = 0; i < block_width * block_height; i += 8) { - i32 j = i/8; - // The order bitplanes are stored in depends on the compressor version - i32 shift_amount; - if (compressor_version == 1) { - shift_amount = (running_bit_index == 0) ? 15 : running_bit_index - 1; // bitplanes are stored sign, lsb ... msb - } else { - shift_amount = 15 - running_bit_index; // bitplanes are stored sign, msb ... lsb - } - u8 b = bitplane[j]; - if (b == 0) continue; + // unpack bitplanes + i32 compressed_bitplane_index = 0; + arena_align(temp_memory.arena, 32); + u16* coeff_buffer = (u16*)arena_push_size(temp_memory.arena, coeff_buffer_size); + memset(coeff_buffer, 0, coeff_buffer_size); + memset(out_buffer, 0, coeff_buffer_size); + + { + u32 running_bit_index = 0; + i32 running_coeff_index = 0; + u32 bitmasks_copy[3]; + memcpy(bitmasks_copy, bitmasks, sizeof(bitmasks)); + for (i32 bitplane_index = 0; bitplane_index < total_mask_bits; ++bitplane_index) { + u8* bitplane = decompressed_buffer + (bitplane_index * bytes_per_bitplane); + + // horribly complicated 'for loop'-style iteration, needed because v1 and v2 store bitplanes in a different order + // v1: iterate over the bitplanes (each coefficient separately) + // v2: for each bitplane, alternate the coefficients (striped) + // TODO: refactor? + if (compressor_version == 1) { + // find next bit for current coeff + for (;;) { + if (running_coeff_index >= coeff_count) { + panic("too many bitplanes"); + } + u16 bitmask = bitmasks_copy[running_coeff_index]; + if (bitmask) { + running_bit_index = bit_scan_forward(bitmask); + ASSERT(running_bit_index < 16); + bitmasks_copy[running_coeff_index] &= ~(1 << running_bit_index); + break; // success + } else { + ++running_coeff_index; + } + } + } else { + // compressor version 2: alternating coefficients + for (;;) { + if (running_bit_index >= 16) { + panic("too many bitplanes"); + } + if (running_coeff_index < coeff_count) { + u16 bitmask = bitmasks_copy[running_coeff_index]; + if (bitmask & (1 << running_bit_index)) { + bitmasks_copy[running_coeff_index] &= ~(1 << running_bit_index); + break; // success + } else { + ++running_coeff_index; // keep looking + } + } else { + running_coeff_index = 0; + ++running_bit_index; // keep looking + } + } + } + // Now we figured out which coeff and bit number this bitplane belongs to + + u16* current_coeff_buffer = coeff_buffer + (running_coeff_index * (block_width * block_height)); + u16* current_out_buffer = (u16*)out_buffer + (running_coeff_index * (block_width * block_height)); + + // Do the bitplane unpacking + for (i32 i = 0; i < block_width * block_height; i += 8) { + i32 j = i/8; + // The order bitplanes are stored in depends on the compressor version + i32 shift_amount; + if (compressor_version == 1) { + shift_amount = (running_bit_index == 0) ? 15 : running_bit_index - 1; // bitplanes are stored sign, lsb ... msb + } else { + shift_amount = 15 - running_bit_index; // bitplanes are stored sign, msb ... lsb + } + u8 b = bitplane[j]; + if (b == 0) continue; #if !defined(__SSE2__) - // Non-SIMD version - current_coeff_buffer[i+0] |= ((b >> 0) & 1) << shift_amount; - current_coeff_buffer[i+1] |= ((b >> 1) & 1) << shift_amount; - current_coeff_buffer[i+2] |= ((b >> 2) & 1) << shift_amount; - current_coeff_buffer[i+3] |= ((b >> 3) & 1) << shift_amount; - current_coeff_buffer[i+4] |= ((b >> 4) & 1) << shift_amount; - current_coeff_buffer[i+5] |= ((b >> 5) & 1) << shift_amount; - current_coeff_buffer[i+6] |= ((b >> 6) & 1) << shift_amount; - current_coeff_buffer[i+7] |= ((b >> 7) & 1) << shift_amount; + // Non-SIMD version + current_coeff_buffer[i+0] |= ((b >> 0) & 1) << shift_amount; + current_coeff_buffer[i+1] |= ((b >> 1) & 1) << shift_amount; + current_coeff_buffer[i+2] |= ((b >> 2) & 1) << shift_amount; + current_coeff_buffer[i+3] |= ((b >> 3) & 1) << shift_amount; + current_coeff_buffer[i+4] |= ((b >> 4) & 1) << shift_amount; + current_coeff_buffer[i+5] |= ((b >> 5) & 1) << shift_amount; + current_coeff_buffer[i+6] |= ((b >> 6) & 1) << shift_amount; + current_coeff_buffer[i+7] |= ((b >> 7) & 1) << shift_amount; #else - // This SIMD implementation is ~20% faster compared to the simple version above. + // This SIMD implementation is ~20% faster compared to the simple version above. // Can it be made faster? __m128i* dst = (__m128i*) (current_coeff_buffer+i); uint64_t t = bswap_64(((0x8040201008040201ULL*b) & 0x8080808080808080ULL) >> 7); @@ -2614,377 +2617,377 @@ bool isyntax_hulsken_decompress(u8* compressed, size_t compressed_size, i32 bloc __m128i result = _mm_or_si128(*dst, masks); *dst = result; #endif - DUMMY_STATEMENT; - } - - // finish iterating: basically, this is the '++i' part of the 'for loop' that got complicated because - // the v1 and v2 compressors store bitplanes in a different order - // TODO: refactor? - if (compressor_version == 2) { - ++running_coeff_index; - } - } - } - - // Reshuffle 4x2 snake-order and convert signed magnitude to twos complement - for (i32 coeff_index = 0; coeff_index < coeff_count; ++coeff_index) { - u16 bitmask = bitmasks[coeff_index]; - u16* current_coeff_buffer = coeff_buffer + (coeff_index * (block_width * block_height)); - u16* current_out_buffer = (u16*)out_buffer + (coeff_index * (block_width * block_height)); - if (bitmask > 0) { - // Reshuffle snake-order - i32 area_stride_x = block_width / 4; - for (i32 area4x4_index = 0; area4x4_index < ((block_width * block_height) / 16); ++area4x4_index) { - i32 area_base_index = area4x4_index * 16; - i32 area_x = (area4x4_index % area_stride_x) * 4; - i32 area_y = (area4x4_index / area_stride_x) * 4; - - u64 area_y0 = *(u64*)¤t_coeff_buffer[area_base_index]; - u64 area_y1 = *(u64*)¤t_coeff_buffer[area_base_index+4]; - u64 area_y2 = *(u64*)¤t_coeff_buffer[area_base_index+8]; - u64 area_y3 = *(u64*)¤t_coeff_buffer[area_base_index+12]; - - *(u64*)(current_out_buffer + (area_y + 0) * block_width + area_x) = area_y0; - *(u64*)(current_out_buffer + (area_y + 1) * block_width + area_x) = area_y1; - *(u64*)(current_out_buffer + (area_y + 2) * block_width + area_x) = area_y2; - *(u64*)(current_out_buffer + (area_y + 3) * block_width + area_x) = area_y3; - } + DUMMY_STATEMENT; + } + + // finish iterating: basically, this is the '++i' part of the 'for loop' that got complicated because + // the v1 and v2 compressors store bitplanes in a different order + // TODO: refactor? + if (compressor_version == 2) { + ++running_coeff_index; + } + } + } - // Convert signed magnitude to twos complement (ex. 0x8002 becomes -2) - signed_magnitude_to_twos_complement_16_block(current_out_buffer, block_width * block_height); - } - } + // Reshuffle 4x2 snake-order and convert signed magnitude to twos complement + for (i32 coeff_index = 0; coeff_index < coeff_count; ++coeff_index) { + u16 bitmask = bitmasks[coeff_index]; + u16* current_coeff_buffer = coeff_buffer + (coeff_index * (block_width * block_height)); + u16* current_out_buffer = (u16*)out_buffer + (coeff_index * (block_width * block_height)); + if (bitmask > 0) { + // Reshuffle snake-order + i32 area_stride_x = block_width / 4; + for (i32 area4x4_index = 0; area4x4_index < ((block_width * block_height) / 16); ++area4x4_index) { + i32 area_base_index = area4x4_index * 16; + i32 area_x = (area4x4_index % area_stride_x) * 4; + i32 area_y = (area4x4_index / area_stride_x) * 4; + + u64 area_y0 = *(u64*)¤t_coeff_buffer[area_base_index]; + u64 area_y1 = *(u64*)¤t_coeff_buffer[area_base_index+4]; + u64 area_y2 = *(u64*)¤t_coeff_buffer[area_base_index+8]; + u64 area_y3 = *(u64*)¤t_coeff_buffer[area_base_index+12]; + + *(u64*)(current_out_buffer + (area_y + 0) * block_width + area_x) = area_y0; + *(u64*)(current_out_buffer + (area_y + 1) * block_width + area_x) = area_y1; + *(u64*)(current_out_buffer + (area_y + 2) * block_width + area_x) = area_y2; + *(u64*)(current_out_buffer + (area_y + 3) * block_width + area_x) = area_y3; + } + + // Convert signed magnitude to twos complement (ex. 0x8002 becomes -2) + signed_magnitude_to_twos_complement_16_block(current_out_buffer, block_width * block_height); + } + } - release_temp_memory(&temp_memory); // frees coeff_buffer and decompressed_buffer - return true; + release_temp_memory(&temp_memory); // frees coeff_buffer and decompressed_buffer + return true; } static inline i32 get_first_valid_coef_pixel(i32 scale) { - i32 result = (PER_LEVEL_PADDING << scale) - (PER_LEVEL_PADDING - 1); - return result; + i32 result = (PER_LEVEL_PADDING << scale) - (PER_LEVEL_PADDING - 1); + return result; } static inline i32 get_first_valid_ll_pixel(i32 scale) { - i32 result = get_first_valid_coef_pixel(scale) + (1 << scale); - return result; + i32 result = get_first_valid_coef_pixel(scale) + (1 << scale); + return result; } i32 isyntax_get_chunk_codeblocks_per_color_for_level(i32 level, bool has_ll) { - i32 rel_level = level % 3; - i32 codeblock_count; - if (rel_level == 0) { - codeblock_count = 1; - } else if (rel_level == 1) { - codeblock_count = 1 + 4; - } else { - codeblock_count = 1 + 4 + 16; - } - if (has_ll) ++codeblock_count; - return codeblock_count; + i32 rel_level = level % 3; + i32 codeblock_count; + if (rel_level == 0) { + codeblock_count = 1; + } else if (rel_level == 1) { + codeblock_count = 1 + 4; + } else { + codeblock_count = 1 + 4 + 16; + } + if (has_ll) ++codeblock_count; + return codeblock_count; } // Dump codeblock info from block header to a .csv file static void isyntax_dump_block_header(isyntax_image_t* wsi_image, const char* filename) { - if (filename == NULL) { - filename = "test_block_header.csv"; - } - FILE* test_block_header_fp = fopen(filename, "wb"); - if (test_block_header_fp) { - fprintf(test_block_header_fp, "x_coordinate,y_coordinate,color_component,scale,coefficient,block_data_offset,block_data_size,block_header_template_id\n"); + if (filename == NULL) { + filename = "test_block_header.csv"; + } + FILE* test_block_header_fp = fopen(filename, "wb"); + if (test_block_header_fp) { + fprintf(test_block_header_fp, "x_coordinate,y_coordinate,color_component,scale,coefficient,block_data_offset,block_data_size,block_header_template_id\n"); - for (i32 i = 0; i < wsi_image->codeblock_count; i += 1/*21*3*/) { - isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; - fprintf(test_block_header_fp, "%d,%d,%d,%d,%d,%lld,%lld,%d\n", + for (i32 i = 0; i < wsi_image->codeblock_count; i += 1/*21*3*/) { + isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; + fprintf(test_block_header_fp, "%d,%d,%d,%d,%d,%lld,%lld,%d\n", // codeblock->x_adjusted, // codeblock->y_adjusted, - codeblock->x_coordinate - wsi_image->offset_x, - codeblock->y_coordinate - wsi_image->offset_y, - codeblock->color_component, - codeblock->scale, - codeblock->coefficient, - codeblock->block_data_offset, - codeblock->block_size, - codeblock->block_header_template_id); + codeblock->x_coordinate - wsi_image->offset_x, + codeblock->y_coordinate - wsi_image->offset_y, + codeblock->color_component, + codeblock->scale, + codeblock->coefficient, + codeblock->block_data_offset, + codeblock->block_size, + codeblock->block_header_template_id); // if (codeblock->scale == wsi_image->level_count-1) { // i+=3; // skip extra LL blocks; // } - } + } - fclose(test_block_header_fp); - } + fclose(test_block_header_fp); + } } // Set the work queue to submit parallel jobs to // TODO(pvalkema): remove this? needs rethinking void isyntax_set_work_queue(isyntax_t* isyntax, work_queue_t* work_queue) { - isyntax->work_submission_queue = work_queue; + isyntax->work_submission_queue = work_queue; } bool isyntax_open(isyntax_t* isyntax, const char* filename, bool init_allocators) { - console_print_verbose("Attempting to open iSyntax: %s\n", filename); - ASSERT(isyntax); - - int ret = 0; (void)ret; - file_stream_t fp = file_stream_open_for_reading(filename); - bool success = false; - - char* read_buffer = NULL; - isyntax_seektable_codeblock_header_t* seektable = NULL; - isyntax_data_chunk_t* data_chunks_memory = NULL; - - if (0) { failed: - if (fp) file_stream_close(fp); - if (read_buffer != NULL) free(read_buffer); - if (seektable != NULL) free(seektable); - isyntax_image_t* wsi_image = isyntax->images + isyntax->wsi_image_index; - if (wsi_image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { - if (wsi_image->data_chunks != NULL) free(wsi_image->data_chunks); - for (i32 i = 0; i < wsi_image->level_count; ++i) { - isyntax_level_t* level = wsi_image->levels + i; - if (level->tiles != NULL) free(level->tiles); - } - } - return success; - } - - if (fp) { - i64 filesize = file_stream_get_filesize(fp); - if (filesize > 0) { - isyntax->filesize = filesize; - - // https://www.openpathology.philips.com/wp-content/uploads/isyntax/4522%20207%2043941_2020_04_24%20Pathology%20iSyntax%20image%20format.pdf - // Layout of an iSyntax file: - // XML Header | End of Table (EOT) marker, 3 bytes "\r\n\x04" | Seektable (optional) | Codeblocks - - // Read the XML header - // We don't know the length of the XML header, so we just read 'enough' data in a chunk and hope that we get it - // (and if not, read some more data until we get it) - - i64 load_begin = get_clock(); - i64 io_begin = get_clock(); - i64 io_ticks_elapsed = 0; - i64 parse_begin = get_clock(); - i64 parse_ticks_elapsed = 0; - - size_t read_size = MEGABYTES(1); - read_buffer = malloc(read_size); - size_t bytes_read = file_stream_read(read_buffer, read_size, fp); - io_ticks_elapsed += (get_clock() - io_begin); - - if (bytes_read < 3) { - goto failed; - } - bool are_there_bytes_left = (bytes_read == read_size); - - // find EOT candidates, 3 bytes "\r\n\x04" - i64 header_length = 0; - i64 isyntax_data_offset = 0; // Offset of either the Seektable, or the Codeblocks segment of the iSyntax file. - i64 bytes_read_from_data_offset_in_last_chunk = 0; + console_print_verbose("Attempting to open iSyntax: %s\n", filename); + ASSERT(isyntax); + + int ret = 0; (void)ret; + file_stream_t fp = file_stream_open_for_reading(filename); + bool success = false; + + char* read_buffer = NULL; + isyntax_seektable_codeblock_header_t* seektable = NULL; + isyntax_data_chunk_t* data_chunks_memory = NULL; + + if (0) { failed: + if (fp) file_stream_close(fp); + if (read_buffer != NULL) free(read_buffer); + if (seektable != NULL) free(seektable); + isyntax_image_t* wsi_image = isyntax->images + isyntax->wsi_image_index; + if (wsi_image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { + if (wsi_image->data_chunks != NULL) free(wsi_image->data_chunks); + for (i32 i = 0; i < wsi_image->level_count; ++i) { + isyntax_level_t* level = wsi_image->levels + i; + if (level->tiles != NULL) free(level->tiles); + } + } + return success; + } - i32 chunk_index = 0; - for (;; ++chunk_index) { + if (fp) { + i64 filesize = file_stream_get_filesize(fp); + if (filesize > 0) { + isyntax->filesize = filesize; + + // https://www.openpathology.philips.com/wp-content/uploads/isyntax/4522%20207%2043941_2020_04_24%20Pathology%20iSyntax%20image%20format.pdf + // Layout of an iSyntax file: + // XML Header | End of Table (EOT) marker, 3 bytes "\r\n\x04" | Seektable (optional) | Codeblocks + + // Read the XML header + // We don't know the length of the XML header, so we just read 'enough' data in a chunk and hope that we get it + // (and if not, read some more data until we get it) + + i64 load_begin = get_clock(); + i64 io_begin = get_clock(); + i64 io_ticks_elapsed = 0; + i64 parse_begin = get_clock(); + i64 parse_ticks_elapsed = 0; + + size_t read_size = MEGABYTES(1); + read_buffer = malloc(read_size); + size_t bytes_read = file_stream_read(read_buffer, read_size, fp); + io_ticks_elapsed += (get_clock() - io_begin); + + if (bytes_read < 3) { + goto failed; + } + bool are_there_bytes_left = (bytes_read == read_size); + + // find EOT candidates, 3 bytes "\r\n\x04" + i64 header_length = 0; + i64 isyntax_data_offset = 0; // Offset of either the Seektable, or the Codeblocks segment of the iSyntax file. + i64 bytes_read_from_data_offset_in_last_chunk = 0; + + i32 chunk_index = 0; + for (;; ++chunk_index) { // console_print_verbose("iSyntax: reading XML header chunk %d\n", chunk_index); - i64 chunk_length = 0; - bool match = false; - char* pos = read_buffer; - i64 offset = 0; - char* marker = (char*)memchr(read_buffer, '\x04', bytes_read); - if (marker) { - offset = marker - read_buffer; - match = true; - chunk_length = offset; - header_length += chunk_length; - isyntax_data_offset = header_length + 1; - i64 data_offset_in_last_chunk = offset + 1; - bytes_read_from_data_offset_in_last_chunk = (i64)bytes_read - data_offset_in_last_chunk; - } - if (match) { - // We found the end of the XML header. This is the last chunk to process. - if (!(header_length > 0 && header_length < isyntax->filesize)) { - goto failed; - } - - parse_begin = get_clock(); - if (!isyntax_parse_xml_header(isyntax, read_buffer, chunk_length, true)) { - goto failed; - } - parse_ticks_elapsed += (get_clock() - parse_begin); + i64 chunk_length = 0; + bool match = false; + char* pos = read_buffer; + i64 offset = 0; + char* marker = (char*)memchr(read_buffer, '\x04', bytes_read); + if (marker) { + offset = marker - read_buffer; + match = true; + chunk_length = offset; + header_length += chunk_length; + isyntax_data_offset = header_length + 1; + i64 data_offset_in_last_chunk = offset + 1; + bytes_read_from_data_offset_in_last_chunk = (i64)bytes_read - data_offset_in_last_chunk; + } + if (match) { + // We found the end of the XML header. This is the last chunk to process. + if (!(header_length > 0 && header_length < isyntax->filesize)) { + goto failed; + } + + parse_begin = get_clock(); + if (!isyntax_parse_xml_header(isyntax, read_buffer, chunk_length, true)) { + goto failed; + } + parse_ticks_elapsed += (get_clock() - parse_begin); // console_print("iSyntax: the XML header is %u bytes, or %g%% of the total file size\n", header_length, (float)((float)header_length * 100.0f) / isyntax->filesize); // console_print(" I/O time: %g seconds\n", get_seconds_elapsed(0, io_ticks_elapsed)); // console_print(" Parsing time: %g seconds\n", get_seconds_elapsed(0, parse_ticks_elapsed)); // console_print(" Total loading time: %g seconds\n", get_seconds_elapsed(load_begin, get_clock())); - break; - } else { - // We didn't find the end of the XML header. We need to read more chunks to find it. - // (Or, we reached the end of the file unexpectedly, which is an error.) - chunk_length = read_size; - header_length += chunk_length; - if (are_there_bytes_left) { - - parse_begin = get_clock(); - if (!isyntax_parse_xml_header(isyntax, read_buffer, chunk_length, false)) { - goto failed; - } - parse_ticks_elapsed += (get_clock() - parse_begin); - - io_begin = get_clock(); - bytes_read = file_stream_read(read_buffer, read_size, fp); // read the next chunk - io_ticks_elapsed += (get_clock() - io_begin); - - are_there_bytes_left = (bytes_read == read_size); - continue; - } else { - console_print_error("iSyntax parsing error: didn't find the end of the XML header (unexpected end of file)\n"); - goto failed; - } - } - } - - if (isyntax->mpp_x <= 0.0f || isyntax->mpp_y <= 0.0f) { - isyntax->mpp_x = 1.0f; // should usually be 0.25; zero or below can never be right - isyntax->mpp_y = 1.0f; - isyntax->is_mpp_known = false; - } - - isyntax->block_width = isyntax->block_header_templates[0].block_width; - isyntax->block_height = isyntax->block_header_templates[0].block_height; - isyntax->tile_width = isyntax->block_width * 2; // tile dimension AFTER inverse wavelet transform - isyntax->tile_height = isyntax->block_height * 2; - - isyntax_image_t* wsi_image = isyntax->images + isyntax->wsi_image_index; - if (wsi_image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { - - i32 block_width = isyntax->block_width; - i32 block_height = isyntax->block_height; - i32 tile_width = isyntax->tile_width; - i32 tile_height = isyntax->tile_height; - - i32 num_levels = wsi_image->level_count; - ASSERT(num_levels >= 1); - i32 grid_width = ((wsi_image->width + (block_width << num_levels) - 1) / (block_width << num_levels)) << (num_levels - 1); - i32 grid_height = ((wsi_image->height + (block_height << num_levels) - 1) / (block_height << num_levels)) << (num_levels - 1); - - u64 h_coeff_tile_count = 0; // number of tiles with LH/HL/HH coefficients - i32 base_level_tile_count = grid_height * grid_width; - for (i32 scale = 0; scale < wsi_image->level_count; ++scale) { - isyntax_level_t* level = wsi_image->levels + scale; - level->tile_count = base_level_tile_count >> (scale * 2); - h_coeff_tile_count += level->tile_count; - level->scale = scale; - level->width_in_tiles = grid_width >> scale; - level->height_in_tiles = grid_height >> scale; - level->width_in_pixels = wsi_image->width >> scale; - level->height_in_pixels = wsi_image->height >> scale; - level->downsample_factor = (float)(1 << scale); - level->um_per_pixel_x = isyntax->mpp_x * level->downsample_factor; - level->um_per_pixel_y = isyntax->mpp_y * level->downsample_factor; - level->x_tile_side_in_um = (float)tile_width * level->um_per_pixel_x; - level->y_tile_side_in_um = (float)tile_height * level->um_per_pixel_y; - } - - // When recursively decoding the tiles, at each iteration the image is slightly offset - // to the top left. - // The shift corresponds to the per level padding added for the wavelet transform: - // ((3 << (scale-1)) - 2) - // Put another way: the highest (zoomed out levels) are shifted the to the bottom right - // (this is also reflected in the x and y coordinates of the codeblocks in the iSyntax header). - // Level 0 has no - for (i32 scale = 1; scale < wsi_image->level_count; ++scale) { - isyntax_level_t* level = wsi_image->levels + scale; - level->origin_offset_in_pixels = get_first_valid_coef_pixel(scale - 1); - float offset_in_um_x = (float)level->origin_offset_in_pixels * wsi_image->levels[0].um_per_pixel_x; - float offset_in_um_y = (float)level->origin_offset_in_pixels * wsi_image->levels[0].um_per_pixel_y; - level->origin_offset = (v2f){offset_in_um_x, offset_in_um_y}; - } - - // The highest level has LL tiles in addition to LH/HL/HH tiles - i64 ll_coeff_tile_count = base_level_tile_count >> ((num_levels - 1) * 2); - i64 total_coeff_tile_count = h_coeff_tile_count + ll_coeff_tile_count; - i64 total_codeblock_count = total_coeff_tile_count * 3; // for 3 color channels - - for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { - isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; - - // Calculate adjusted codeblock coordinates so that they fit the origin of the image - codeblock->x_adjusted = (i32)codeblock->x_coordinate - wsi_image->offset_x; - codeblock->y_adjusted = (i32)codeblock->y_coordinate - wsi_image->offset_y; - - // Calculate the block ID (= index into the seektable) - // adapted from extract_block_header.py - bool is_ll = codeblock->coefficient == 0; - u32 block_id = 0; - i32 maxscale = is_ll ? codeblock->scale + 1 : codeblock->scale; - for (i32 scale = 0; scale < maxscale; ++scale) { - block_id += wsi_image->levels[scale].tile_count; - } - - i32 offset; - if (is_ll) { - offset = get_first_valid_ll_pixel(codeblock->scale); - } else { - offset = get_first_valid_coef_pixel(codeblock->scale); - } - i32 x = codeblock->x_adjusted - offset; - i32 y = codeblock->y_adjusted - offset; + break; + } else { + // We didn't find the end of the XML header. We need to read more chunks to find it. + // (Or, we reached the end of the file unexpectedly, which is an error.) + chunk_length = read_size; + header_length += chunk_length; + if (are_there_bytes_left) { + + parse_begin = get_clock(); + if (!isyntax_parse_xml_header(isyntax, read_buffer, chunk_length, false)) { + goto failed; + } + parse_ticks_elapsed += (get_clock() - parse_begin); + + io_begin = get_clock(); + bytes_read = file_stream_read(read_buffer, read_size, fp); // read the next chunk + io_ticks_elapsed += (get_clock() - io_begin); + + are_there_bytes_left = (bytes_read == read_size); + continue; + } else { + console_print_error("iSyntax parsing error: didn't find the end of the XML header (unexpected end of file)\n"); + goto failed; + } + } + } + + if (isyntax->mpp_x <= 0.0f || isyntax->mpp_y <= 0.0f) { + isyntax->mpp_x = 1.0f; // should usually be 0.25; zero or below can never be right + isyntax->mpp_y = 1.0f; + isyntax->is_mpp_known = false; + } + + isyntax->block_width = isyntax->block_header_templates[0].block_width; + isyntax->block_height = isyntax->block_header_templates[0].block_height; + isyntax->tile_width = isyntax->block_width * 2; // tile dimension AFTER inverse wavelet transform + isyntax->tile_height = isyntax->block_height * 2; + + isyntax_image_t* wsi_image = isyntax->images + isyntax->wsi_image_index; + if (wsi_image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { + + i32 block_width = isyntax->block_width; + i32 block_height = isyntax->block_height; + i32 tile_width = isyntax->tile_width; + i32 tile_height = isyntax->tile_height; + + i32 num_levels = wsi_image->level_count; + ASSERT(num_levels >= 1); + i32 grid_width = ((wsi_image->width + (block_width << num_levels) - 1) / (block_width << num_levels)) << (num_levels - 1); + i32 grid_height = ((wsi_image->height + (block_height << num_levels) - 1) / (block_height << num_levels)) << (num_levels - 1); + + u64 h_coeff_tile_count = 0; // number of tiles with LH/HL/HH coefficients + i32 base_level_tile_count = grid_height * grid_width; + for (i32 scale = 0; scale < wsi_image->level_count; ++scale) { + isyntax_level_t* level = wsi_image->levels + scale; + level->tile_count = base_level_tile_count >> (scale * 2); + h_coeff_tile_count += level->tile_count; + level->scale = scale; + level->width_in_tiles = grid_width >> scale; + level->height_in_tiles = grid_height >> scale; + level->width_in_pixels = wsi_image->width >> scale; + level->height_in_pixels = wsi_image->height >> scale; + level->downsample_factor = (float)(1 << scale); + level->um_per_pixel_x = isyntax->mpp_x * level->downsample_factor; + level->um_per_pixel_y = isyntax->mpp_y * level->downsample_factor; + level->x_tile_side_in_um = (float)tile_width * level->um_per_pixel_x; + level->y_tile_side_in_um = (float)tile_height * level->um_per_pixel_y; + } + + // When recursively decoding the tiles, at each iteration the image is slightly offset + // to the top left. + // The shift corresponds to the per level padding added for the wavelet transform: + // ((3 << (scale-1)) - 2) + // Put another way: the highest (zoomed out levels) are shifted the to the bottom right + // (this is also reflected in the x and y coordinates of the codeblocks in the iSyntax header). + // Level 0 has no + for (i32 scale = 1; scale < wsi_image->level_count; ++scale) { + isyntax_level_t* level = wsi_image->levels + scale; + level->origin_offset_in_pixels = get_first_valid_coef_pixel(scale - 1); + float offset_in_um_x = (float)level->origin_offset_in_pixels * wsi_image->levels[0].um_per_pixel_x; + float offset_in_um_y = (float)level->origin_offset_in_pixels * wsi_image->levels[0].um_per_pixel_y; + level->origin_offset = (v2f){offset_in_um_x, offset_in_um_y}; + } + + // The highest level has LL tiles in addition to LH/HL/HH tiles + i64 ll_coeff_tile_count = base_level_tile_count >> ((num_levels - 1) * 2); + i64 total_coeff_tile_count = h_coeff_tile_count + ll_coeff_tile_count; + i64 total_codeblock_count = total_coeff_tile_count * 3; // for 3 color channels + + for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { + isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; + + // Calculate adjusted codeblock coordinates so that they fit the origin of the image + codeblock->x_adjusted = (i32)codeblock->x_coordinate - wsi_image->offset_x; + codeblock->y_adjusted = (i32)codeblock->y_coordinate - wsi_image->offset_y; + + // Calculate the block ID (= index into the seektable) + // adapted from extract_block_header.py + bool is_ll = codeblock->coefficient == 0; + u32 block_id = 0; + i32 maxscale = is_ll ? codeblock->scale + 1 : codeblock->scale; + for (i32 scale = 0; scale < maxscale; ++scale) { + block_id += wsi_image->levels[scale].tile_count; + } + + i32 offset; + if (is_ll) { + offset = get_first_valid_ll_pixel(codeblock->scale); + } else { + offset = get_first_valid_coef_pixel(codeblock->scale); + } + i32 x = codeblock->x_adjusted - offset; + i32 y = codeblock->y_adjusted - offset; // codeblock->x_adjusted = x; // codeblock->y_adjusted = y; - codeblock->block_x = x / (tile_width << codeblock->scale); - codeblock->block_y = y / (tile_height << codeblock->scale); - - i32 grid_stride = grid_width >> codeblock->scale; - block_id += codeblock->block_y * grid_stride + codeblock->block_x; - - i32 tiles_per_color = total_coeff_tile_count; - block_id += codeblock->color_component * tiles_per_color; - codeblock->block_id = block_id; - } - - io_begin = get_clock(); // for performance measurement - file_stream_set_pos(fp, isyntax_data_offset); - if (wsi_image->header_codeblocks_are_partial) { - // The seektable is required to be present, because the block header table did not contain all information. - isyntax_dicom_tag_header_t seektable_header_tag = {0}; - file_stream_read(&seektable_header_tag, sizeof(isyntax_dicom_tag_header_t), fp); - - io_ticks_elapsed += (get_clock() - io_begin); - parse_begin = get_clock(); - - if (seektable_header_tag.group == 0x301D && seektable_header_tag.element == 0x2015) { - i32 seektable_size = seektable_header_tag.size; - if (seektable_size < 0) { - // We need to guess the size... - ASSERT(wsi_image->codeblock_count > 0); - seektable_size = sizeof(isyntax_seektable_codeblock_header_t) * wsi_image->codeblock_count; - } - seektable = (isyntax_seektable_codeblock_header_t*) malloc(seektable_size); - file_stream_read(seektable, seektable_size, fp); - - // Now fill in the missing data. - // NOTE: The number of codeblock entries in the seektable is much greater than the number of - // codeblocks that *actually* exist in the file. This means that we have to discard many of - // the seektable entries. - // Luckily, we can easily identify the entries that need to be discarded. - // (They have the data offset (and data size) set to 0.) - i32 seektable_entry_count = seektable_size / sizeof(isyntax_seektable_codeblock_header_t); - - for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { - isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; - if (codeblock->block_id > seektable_entry_count) { - ASSERT(!"block ID out of bounds"); - goto failed; - } - isyntax_seektable_codeblock_header_t* seektable_entry = seektable + codeblock->block_id; - ASSERT(seektable_entry->block_data_offset_header.group == 0x301D); - ASSERT(seektable_entry->block_data_offset_header.element == 0x2010); - codeblock->block_data_offset = seektable_entry->block_data_offset; - codeblock->block_size = seektable_entry->block_size; + codeblock->block_x = x / (tile_width << codeblock->scale); + codeblock->block_y = y / (tile_height << codeblock->scale); + + i32 grid_stride = grid_width >> codeblock->scale; + block_id += codeblock->block_y * grid_stride + codeblock->block_x; + + i32 tiles_per_color = total_coeff_tile_count; + block_id += codeblock->color_component * tiles_per_color; + codeblock->block_id = block_id; + } + + io_begin = get_clock(); // for performance measurement + file_stream_set_pos(fp, isyntax_data_offset); + if (wsi_image->header_codeblocks_are_partial) { + // The seektable is required to be present, because the block header table did not contain all information. + isyntax_dicom_tag_header_t seektable_header_tag = {0}; + file_stream_read(&seektable_header_tag, sizeof(isyntax_dicom_tag_header_t), fp); + + io_ticks_elapsed += (get_clock() - io_begin); + parse_begin = get_clock(); + + if (seektable_header_tag.group == 0x301D && seektable_header_tag.element == 0x2015) { + i32 seektable_size = seektable_header_tag.size; + if (seektable_size < 0) { + // We need to guess the size... + ASSERT(wsi_image->codeblock_count > 0); + seektable_size = sizeof(isyntax_seektable_codeblock_header_t) * wsi_image->codeblock_count; + } + seektable = (isyntax_seektable_codeblock_header_t*) malloc(seektable_size); + file_stream_read(seektable, seektable_size, fp); + + // Now fill in the missing data. + // NOTE: The number of codeblock entries in the seektable is much greater than the number of + // codeblocks that *actually* exist in the file. This means that we have to discard many of + // the seektable entries. + // Luckily, we can easily identify the entries that need to be discarded. + // (They have the data offset (and data size) set to 0.) + i32 seektable_entry_count = seektable_size / sizeof(isyntax_seektable_codeblock_header_t); + + for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { + isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; + if (codeblock->block_id > seektable_entry_count) { + ASSERT(!"block ID out of bounds"); + goto failed; + } + isyntax_seektable_codeblock_header_t* seektable_entry = seektable + codeblock->block_id; + ASSERT(seektable_entry->block_data_offset_header.group == 0x301D); + ASSERT(seektable_entry->block_data_offset_header.element == 0x2010); + codeblock->block_data_offset = seektable_entry->block_data_offset; + codeblock->block_size = seektable_entry->block_size; #if 0 - // Debug test: + // Debug test: // Decompress codeblocks in the seektable. if (1 || i == wsi_image->codeblock_count_per_color-1) { // Only parse 'non-empty'/'background' codeblocks @@ -3012,221 +3015,221 @@ bool isyntax_open(isyntax_t* isyntax, const char* filename, bool init_allocators } #endif - } - free(seektable); - seektable = NULL; + } + free(seektable); + seektable = NULL; // isyntax_dump_block_header(wsi_image, "test_block_header.csv"); - // Allocate enough space for the maximum number of codeblock 'chunks' we can expect - // (the actual number of chunks may be lower, because some tiles might not exist) - i32 max_possible_chunk_count = 0; - for (i32 scale = 0; scale <= wsi_image->max_scale; ++scale) { - if ((scale + 1) % 3 == 0 || scale == wsi_image->max_scale) { - isyntax_level_t* level = wsi_image->levels + scale; - max_possible_chunk_count += level->tile_count; - } - } - wsi_image->data_chunks = (isyntax_data_chunk_t*) calloc(1, max_possible_chunk_count * sizeof(isyntax_data_chunk_t)); - - // Create tables for spatial lookup of codeblocks and codeblock chunks from tile coordinates - for (i32 i = 0; i < wsi_image->level_count; ++i) { - isyntax_level_t* level = wsi_image->levels + i; - // NOTE: Tile entry with codeblock_index == 0 will mean there is no codeblock for this tile (empty/background) - level->tiles = (isyntax_tile_t*) calloc(1, level->tile_count * sizeof(isyntax_tile_t)); - } - i32 current_chunk_codeblock_index = 0; - i32 next_chunk_codeblock_index = 0; - i32 current_data_chunk_index = 0; - i32 next_data_chunk_index = 0; - for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { - isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; - if (codeblock->color_component != 0) { - // Don't let color channels 1 and 2 overwrite what was already set - i = next_chunk_codeblock_index; // skip ahead - codeblock = wsi_image->codeblocks + i; - if (i >= wsi_image->codeblock_count) break; - } - // Keep track of where we are in the 'chunk' of codeblocks - if (i == next_chunk_codeblock_index) { - // This codeblock is the top of a new chunk - i32 chunk_codeblock_count_per_color; - if (codeblock->scale == wsi_image->max_scale) { - chunk_codeblock_count_per_color = isyntax_get_chunk_codeblocks_per_color_for_level(codeblock->scale, true); - } else { - chunk_codeblock_count_per_color = 21; - } - current_chunk_codeblock_index = i; - next_chunk_codeblock_index = i + (chunk_codeblock_count_per_color * 3); - current_data_chunk_index = next_data_chunk_index; - if (current_data_chunk_index >= max_possible_chunk_count) { - console_print_error("iSyntax: encountered too many data chunks\n"); - panic(); - } - - isyntax_data_chunk_t* chunk = wsi_image->data_chunks + current_data_chunk_index; - chunk->offset = codeblock->block_data_offset; - // TODO: record cluster size here? - chunk->top_codeblock_index = current_chunk_codeblock_index; - chunk->codeblock_count_per_color = chunk_codeblock_count_per_color; - chunk->scale = codeblock->scale; - ++wsi_image->data_chunk_count; - ++next_data_chunk_index; - } - isyntax_level_t* level = wsi_image->levels + codeblock->scale; - i32 tile_index = codeblock->block_y * level->width_in_tiles + codeblock->block_x; - ASSERT(tile_index < level->tile_count); - level->tiles[tile_index].exists = true; - level->tiles[tile_index].codeblock_index = i; - level->tiles[tile_index].codeblock_chunk_index = current_chunk_codeblock_index; - level->tiles[tile_index].data_chunk_index = current_data_chunk_index; - - } - - parse_ticks_elapsed += (get_clock() - parse_begin); + // Allocate enough space for the maximum number of codeblock 'chunks' we can expect + // (the actual number of chunks may be lower, because some tiles might not exist) + i32 max_possible_chunk_count = 0; + for (i32 scale = 0; scale <= wsi_image->max_scale; ++scale) { + if ((scale + 1) % 3 == 0 || scale == wsi_image->max_scale) { + isyntax_level_t* level = wsi_image->levels + scale; + max_possible_chunk_count += level->tile_count; + } + } + wsi_image->data_chunks = (isyntax_data_chunk_t*) calloc(1, max_possible_chunk_count * sizeof(isyntax_data_chunk_t)); + + // Create tables for spatial lookup of codeblocks and codeblock chunks from tile coordinates + for (i32 i = 0; i < wsi_image->level_count; ++i) { + isyntax_level_t* level = wsi_image->levels + i; + // NOTE: Tile entry with codeblock_index == 0 will mean there is no codeblock for this tile (empty/background) + level->tiles = (isyntax_tile_t*) calloc(1, level->tile_count * sizeof(isyntax_tile_t)); + } + i32 current_chunk_codeblock_index = 0; + i32 next_chunk_codeblock_index = 0; + i32 current_data_chunk_index = 0; + i32 next_data_chunk_index = 0; + for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { + isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; + if (codeblock->color_component != 0) { + // Don't let color channels 1 and 2 overwrite what was already set + i = next_chunk_codeblock_index; // skip ahead + codeblock = wsi_image->codeblocks + i; + if (i >= wsi_image->codeblock_count) break; + } + // Keep track of where we are in the 'chunk' of codeblocks + if (i == next_chunk_codeblock_index) { + // This codeblock is the top of a new chunk + i32 chunk_codeblock_count_per_color; + if (codeblock->scale == wsi_image->max_scale) { + chunk_codeblock_count_per_color = isyntax_get_chunk_codeblocks_per_color_for_level(codeblock->scale, true); + } else { + chunk_codeblock_count_per_color = 21; + } + current_chunk_codeblock_index = i; + next_chunk_codeblock_index = i + (chunk_codeblock_count_per_color * 3); + current_data_chunk_index = next_data_chunk_index; + if (current_data_chunk_index >= max_possible_chunk_count) { + console_print_error("iSyntax: encountered too many data chunks\n"); + panic(); + } + + isyntax_data_chunk_t* chunk = wsi_image->data_chunks + current_data_chunk_index; + chunk->offset = codeblock->block_data_offset; + // TODO: record cluster size here? + chunk->top_codeblock_index = current_chunk_codeblock_index; + chunk->codeblock_count_per_color = chunk_codeblock_count_per_color; + chunk->scale = codeblock->scale; + ++wsi_image->data_chunk_count; + ++next_data_chunk_index; + } + isyntax_level_t* level = wsi_image->levels + codeblock->scale; + i32 tile_index = codeblock->block_y * level->width_in_tiles + codeblock->block_x; + ASSERT(tile_index < level->tile_count); + level->tiles[tile_index].exists = true; + level->tiles[tile_index].codeblock_index = i; + level->tiles[tile_index].codeblock_chunk_index = current_chunk_codeblock_index; + level->tiles[tile_index].data_chunk_index = current_data_chunk_index; + + } + + parse_ticks_elapsed += (get_clock() - parse_begin); // console_print("iSyntax: the seektable is %u bytes, or %g%% of the total file size\n", seektable_size, (float)((float)seektable_size * 100.0f) / isyntax->filesize); // console_print(" I/O time: %g seconds\n", get_seconds_elapsed(0, io_ticks_elapsed)); // console_print(" Parsing time: %g seconds\n", get_seconds_elapsed(0, parse_ticks_elapsed)); - isyntax->loading_time = get_seconds_elapsed(load_begin, get_clock()); - } else { - // seektable invalid - goto failed; - } - } else if (isyntax->data_model_major_version >= 100) { - // Create tables for spatial lookup of codeblocks and codeblock chunks from tile coordinates - for (i32 i = 0; i < wsi_image->level_count; ++i) { - isyntax_level_t* level = wsi_image->levels + i; - // NOTE: Tile entry with codeblock_index == 0 will mean there is no codeblock for this tile (empty/background) - level->tiles = (isyntax_tile_t*) calloc(1, level->tile_count * sizeof(isyntax_tile_t)); - } - - // TODO: refactor code duplication: this is (probably) mostly wrong! - i32 current_chunk_codeblock_index = 0; - i32 next_chunk_codeblock_index = 0; - i32 current_data_chunk_index = 0; - i32 next_data_chunk_index = 0; - for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { - isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; - if (codeblock->color_component != 0) { - // Don't let color channels 1 and 2 overwrite what was already set - i = next_chunk_codeblock_index; // skip ahead - codeblock = wsi_image->codeblocks + i; - if (i >= wsi_image->codeblock_count) break; - } - // Keep track of where we are in the 'chunk' of codeblocks - if (i == next_chunk_codeblock_index) { - // This codeblock is the top of a new chunk - i32 chunk_codeblock_count_per_color; - if (codeblock->scale == wsi_image->max_scale) { - chunk_codeblock_count_per_color = isyntax_get_chunk_codeblocks_per_color_for_level(codeblock->scale, true); - } else { - chunk_codeblock_count_per_color = 21; - } - current_chunk_codeblock_index = i; - next_chunk_codeblock_index = i + (chunk_codeblock_count_per_color * 3); - current_data_chunk_index = next_data_chunk_index; - if (current_data_chunk_index >= wsi_image->data_chunk_count) { - console_print_error("iSyntax: encountered too many data chunks\n"); - panic(); - } - - if (isyntax->data_model_major_version < 100) { - isyntax_data_chunk_t* chunk = wsi_image->data_chunks + current_data_chunk_index; - chunk->offset = codeblock->block_data_offset; - // TODO: record cluster size here - chunk->top_codeblock_index = current_chunk_codeblock_index; - chunk->codeblock_count_per_color = chunk_codeblock_count_per_color; - chunk->scale = codeblock->scale; - ++wsi_image->data_chunk_count; - } - - ++next_data_chunk_index; - } - isyntax_level_t* level = wsi_image->levels + codeblock->scale; - i32 tile_index = codeblock->block_y * level->width_in_tiles + codeblock->block_x; - ASSERT(tile_index < level->tile_count); - level->tiles[tile_index].exists = true; - level->tiles[tile_index].codeblock_index = i; - level->tiles[tile_index].codeblock_chunk_index = current_chunk_codeblock_index; - level->tiles[tile_index].data_chunk_index = current_data_chunk_index; - - } - - parse_ticks_elapsed += (get_clock() - parse_begin); + isyntax->loading_time = get_seconds_elapsed(load_begin, get_clock()); + } else { + // seektable invalid + goto failed; + } + } else if (isyntax->data_model_major_version >= 100) { + // Create tables for spatial lookup of codeblocks and codeblock chunks from tile coordinates + for (i32 i = 0; i < wsi_image->level_count; ++i) { + isyntax_level_t* level = wsi_image->levels + i; + // NOTE: Tile entry with codeblock_index == 0 will mean there is no codeblock for this tile (empty/background) + level->tiles = (isyntax_tile_t*) calloc(1, level->tile_count * sizeof(isyntax_tile_t)); + } + + // TODO: refactor code duplication: this is (probably) mostly wrong! + i32 current_chunk_codeblock_index = 0; + i32 next_chunk_codeblock_index = 0; + i32 current_data_chunk_index = 0; + i32 next_data_chunk_index = 0; + for (i32 i = 0; i < wsi_image->codeblock_count; ++i) { + isyntax_codeblock_t* codeblock = wsi_image->codeblocks + i; + if (codeblock->color_component != 0) { + // Don't let color channels 1 and 2 overwrite what was already set + i = next_chunk_codeblock_index; // skip ahead + codeblock = wsi_image->codeblocks + i; + if (i >= wsi_image->codeblock_count) break; + } + // Keep track of where we are in the 'chunk' of codeblocks + if (i == next_chunk_codeblock_index) { + // This codeblock is the top of a new chunk + i32 chunk_codeblock_count_per_color; + if (codeblock->scale == wsi_image->max_scale) { + chunk_codeblock_count_per_color = isyntax_get_chunk_codeblocks_per_color_for_level(codeblock->scale, true); + } else { + chunk_codeblock_count_per_color = 21; + } + current_chunk_codeblock_index = i; + next_chunk_codeblock_index = i + (chunk_codeblock_count_per_color * 3); + current_data_chunk_index = next_data_chunk_index; + if (current_data_chunk_index >= wsi_image->data_chunk_count) { + console_print_error("iSyntax: encountered too many data chunks\n"); + panic(); + } + + if (isyntax->data_model_major_version < 100) { + isyntax_data_chunk_t* chunk = wsi_image->data_chunks + current_data_chunk_index; + chunk->offset = codeblock->block_data_offset; + // TODO: record cluster size here + chunk->top_codeblock_index = current_chunk_codeblock_index; + chunk->codeblock_count_per_color = chunk_codeblock_count_per_color; + chunk->scale = codeblock->scale; + ++wsi_image->data_chunk_count; + } + + ++next_data_chunk_index; + } + isyntax_level_t* level = wsi_image->levels + codeblock->scale; + i32 tile_index = codeblock->block_y * level->width_in_tiles + codeblock->block_x; + ASSERT(tile_index < level->tile_count); + level->tiles[tile_index].exists = true; + level->tiles[tile_index].codeblock_index = i; + level->tiles[tile_index].codeblock_chunk_index = current_chunk_codeblock_index; + level->tiles[tile_index].data_chunk_index = current_data_chunk_index; + + } + + parse_ticks_elapsed += (get_clock() - parse_begin); // console_print("iSyntax: the seektable is %u bytes, or %g%% of the total file size\n", seektable_size, (float)((float)seektable_size * 100.0f) / isyntax->filesize); // console_print(" I/O time: %g seconds\n", get_seconds_elapsed(0, io_ticks_elapsed)); // console_print(" Parsing time: %g seconds\n", get_seconds_elapsed(0, parse_ticks_elapsed)); - isyntax->loading_time = get_seconds_elapsed(load_begin, get_clock()); - } else { - // non-partial header blocks are not supported - goto failed; - }; - } else { - // non-WSI images are not supported - goto failed; - } - - size_t ll_coeff_block_size = isyntax->block_width * isyntax->block_height * sizeof(icoeff_t); - size_t block_allocator_maximum_capacity_in_blocks = GIGABYTES(32) / ll_coeff_block_size; - size_t ll_coeff_block_allocator_capacity_in_blocks = block_allocator_maximum_capacity_in_blocks / 4; - size_t h_coeff_block_size = ll_coeff_block_size * 3; - size_t h_coeff_block_allocator_capacity_in_blocks = ll_coeff_block_allocator_capacity_in_blocks * 3; - if (init_allocators) { - isyntax->ll_coeff_block_allocator = malloc(sizeof(block_allocator_t)); - isyntax->h_coeff_block_allocator = malloc(sizeof(block_allocator_t)); - *isyntax->ll_coeff_block_allocator = block_allocator_create(ll_coeff_block_size, ll_coeff_block_allocator_capacity_in_blocks, MEGABYTES(256)); - *isyntax->h_coeff_block_allocator = block_allocator_create(h_coeff_block_size, h_coeff_block_allocator_capacity_in_blocks, MEGABYTES(256)); - isyntax->is_block_allocator_owned = true; - } else { - // The caller must inject the allocators after return of isyntax_open(). - isyntax->ll_coeff_block_allocator = NULL; - isyntax->h_coeff_block_allocator = NULL; - isyntax->is_block_allocator_owned = false; - } - - success = true; - - free(read_buffer); - - // Populate debug info. - for (int scale = 0; scale < wsi_image->level_count; ++scale) { - isyntax_level_t* level = &wsi_image->levels[scale]; - for (int tile_y = 0; tile_y < level->height_in_tiles; ++tile_y) { - for (int tile_x = 0; tile_x < level->width_in_tiles; ++tile_x) { - isyntax_tile_t* tile = &level->tiles[level->width_in_tiles * tile_y + tile_x]; - tile->tile_scale = scale; - tile->tile_x = tile_x; - tile->tile_y = tile_y; - } - } - } - } - file_stream_close(fp); - - if (success) { - isyntax->file_handle = open_file_handle_for_simultaneous_access(filename); - if (!isyntax->file_handle) { - console_print_error("Error: Could not reopen file for asynchronous I/O\n"); - success = false; - } - } - } - return success; + isyntax->loading_time = get_seconds_elapsed(load_begin, get_clock()); + } else { + // non-partial header blocks are not supported + goto failed; + }; + } else { + // non-WSI images are not supported + goto failed; + } + + size_t ll_coeff_block_size = isyntax->block_width * isyntax->block_height * sizeof(icoeff_t); + size_t block_allocator_maximum_capacity_in_blocks = GIGABYTES(32) / ll_coeff_block_size; + size_t ll_coeff_block_allocator_capacity_in_blocks = block_allocator_maximum_capacity_in_blocks / 4; + size_t h_coeff_block_size = ll_coeff_block_size * 3; + size_t h_coeff_block_allocator_capacity_in_blocks = ll_coeff_block_allocator_capacity_in_blocks * 3; + if (init_allocators) { + isyntax->ll_coeff_block_allocator = malloc(sizeof(block_allocator_t)); + isyntax->h_coeff_block_allocator = malloc(sizeof(block_allocator_t)); + *isyntax->ll_coeff_block_allocator = block_allocator_create(ll_coeff_block_size, ll_coeff_block_allocator_capacity_in_blocks, MEGABYTES(256)); + *isyntax->h_coeff_block_allocator = block_allocator_create(h_coeff_block_size, h_coeff_block_allocator_capacity_in_blocks, MEGABYTES(256)); + isyntax->is_block_allocator_owned = true; + } else { + // The caller must inject the allocators after return of isyntax_open(). + isyntax->ll_coeff_block_allocator = NULL; + isyntax->h_coeff_block_allocator = NULL; + isyntax->is_block_allocator_owned = false; + } + + success = true; + + free(read_buffer); + + // Populate debug info. + for (int scale = 0; scale < wsi_image->level_count; ++scale) { + isyntax_level_t* level = &wsi_image->levels[scale]; + for (int tile_y = 0; tile_y < level->height_in_tiles; ++tile_y) { + for (int tile_x = 0; tile_x < level->width_in_tiles; ++tile_x) { + isyntax_tile_t* tile = &level->tiles[level->width_in_tiles * tile_y + tile_x]; + tile->tile_scale = scale; + tile->tile_x = tile_x; + tile->tile_y = tile_y; + } + } + } + } + file_stream_close(fp); + + if (success) { + isyntax->file_handle = open_file_handle_for_simultaneous_access(filename); + if (!isyntax->file_handle) { + console_print_error("Error: Could not reopen file for asynchronous I/O\n"); + success = false; + } + } + } + return success; } void isyntax_destroy(isyntax_t* isyntax) { - while (isyntax->refcount > 0) { - platform_sleep(1); - if (isyntax->work_submission_queue) { - do_worker_work(isyntax->work_submission_queue, 0); - } else { - static bool already_printed = false; - if (!already_printed) { - console_print_error("isyntax_destroy(): work_submission_queue not set; refcount = %d, waiting to reach 0\n", isyntax->refcount); - already_printed = true; - } - } - } + while (isyntax->refcount > 0) { + platform_sleep(1); + if (isyntax->work_submission_queue) { + do_worker_work(isyntax->work_submission_queue, 0); + } else { + static bool already_printed = false; + if (!already_printed) { + console_print_error("isyntax_destroy(): work_submission_queue not set; refcount = %d, waiting to reach 0\n", isyntax->refcount); + already_printed = true; + } + } + } if (isyntax->is_block_allocator_owned) { if (isyntax->ll_coeff_block_allocator->is_valid) { block_allocator_destroy(isyntax->ll_coeff_block_allocator); @@ -3235,40 +3238,40 @@ void isyntax_destroy(isyntax_t* isyntax) { block_allocator_destroy(isyntax->h_coeff_block_allocator); } } - if (isyntax->black_dummy_coeff) { - free(isyntax->black_dummy_coeff); - isyntax->black_dummy_coeff = NULL; - } - if (isyntax->white_dummy_coeff) { - free(isyntax->white_dummy_coeff); - isyntax->white_dummy_coeff = NULL; - } - for (i32 image_index = 0; image_index < isyntax->image_count; ++image_index) { - isyntax_image_t* image = isyntax->images + image_index; - if (image->pixels) { - free(image->pixels); - image->pixels = NULL; - } - if (image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { - if (image->codeblocks) { - free(image->codeblocks); - image->codeblocks = NULL; - } - if (image->data_chunks) { - for (i32 i = 0; i < image->data_chunk_count; ++i) { - isyntax_data_chunk_t* chunk = image->data_chunks + i; - if (chunk->data) { - free(chunk->data); - } - } - free(image->data_chunks); - image->data_chunks = NULL; - } - for (i32 i = 0; i < image->level_count; ++i) { - isyntax_level_t* level = image->levels + i; - if (level->tiles) { + if (isyntax->black_dummy_coeff) { + free(isyntax->black_dummy_coeff); + isyntax->black_dummy_coeff = NULL; + } + if (isyntax->white_dummy_coeff) { + free(isyntax->white_dummy_coeff); + isyntax->white_dummy_coeff = NULL; + } + for (i32 image_index = 0; image_index < isyntax->image_count; ++image_index) { + isyntax_image_t* image = isyntax->images + image_index; + if (image->pixels) { + free(image->pixels); + image->pixels = NULL; + } + if (image->image_type == ISYNTAX_IMAGE_TYPE_WSI) { + if (image->codeblocks) { + free(image->codeblocks); + image->codeblocks = NULL; + } + if (image->data_chunks) { + for (i32 i = 0; i < image->data_chunk_count; ++i) { + isyntax_data_chunk_t* chunk = image->data_chunks + i; + if (chunk->data) { + free(chunk->data); + } + } + free(image->data_chunks); + image->data_chunks = NULL; + } + for (i32 i = 0; i < image->level_count; ++i) { + isyntax_level_t* level = image->levels + i; + if (level->tiles) { #if 0 - for (i32 j = 0; j < level->tile_count; ++j) { + for (i32 j = 0; j < level->tile_count; ++j) { isyntax_tile_t* tile = level->tiles + j; for (i32 color = 0; color < 3; ++color) { isyntax_tile_channel_t* channel = tile->color_channels + color; @@ -3277,13 +3280,12 @@ void isyntax_destroy(isyntax_t* isyntax) { } } #endif - free(level->tiles); - level->tiles = NULL; - } - } - } - } - file_handle_close(isyntax->file_handle); + free(level->tiles); + level->tiles = NULL; + } + } + } + } + file_handle_close(isyntax->file_handle); } - diff --git a/src/isyntax/isyntax_reader.c b/src/isyntax/isyntax_reader.c index f2636be..2a5c596 100644 --- a/src/isyntax/isyntax_reader.c +++ b/src/isyntax/isyntax_reader.c @@ -117,7 +117,8 @@ static void isyntax_openslide_load_tile_coefficients_ll_or_h(isyntax_cache_t* ca tile->color_channels[color].coeff_h = (icoeff_t *) block_alloc(&cache->h_coeff_block_allocator); } // TODO(avirodov): fancy allocators, for multiple sequential blocks (aka chunk). Or let OS do the caching. - u8* codeblock_data = malloc(codeblock->block_size); + // Adding 7 safety bytes so bitstream_lsb_read() won't access out of bounds in isyntax_hulsken_decompress(). + u8* codeblock_data = malloc(codeblock->block_size + 7); size_t bytes_read = file_handle_read_at_offset(codeblock_data, isyntax->file_handle, codeblock->block_data_offset, codeblock->block_size); if (!(bytes_read > 0)) { diff --git a/src/isyntax_example.c b/src/isyntax_example.c index 5bba620..01ac9cd 100644 --- a/src/isyntax_example.c +++ b/src/isyntax_example.c @@ -87,4 +87,4 @@ int main(int argc, char** argv) { libisyntax_close(isyntax); return 0; -} +} \ No newline at end of file diff --git a/src/isyntax_to_tiff.c b/src/isyntax_to_tiff.c new file mode 100644 index 0000000..a9644a5 --- /dev/null +++ b/src/isyntax_to_tiff.c @@ -0,0 +1,422 @@ +/* + BSD 2-Clause License + + Copyright (c) 2019-2023, Pieter Valkema, Alexandr Virodov, Jonas Teuwen + + Redistribution and use in source and binary forms, with or without + modification, are permitted provided that the following conditions are met: + + 1. Redistributions of source code must retain the above copyright notice, this + list of conditions and the following disclaimer. + + 2. Redistributions in binary form must reproduce the above copyright notice, + this list of conditions and the following disclaimer in the documentation + and/or other materials provided with the distribution. + + THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" + AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE + IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE + DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE + FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL + DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR + SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER + CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, + OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE + OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +*/ + +#include "libisyntax.h" +#include "tiffio.h" +#include +#include +#include +#include +#include + + +// https://stackoverflow.com/questions/11228855/header-files-for-x86-simd-intrinsics +#if defined(_MSC_VER) +/* Microsoft C/C++-compatible compiler */ + #include +#elif defined(__GNUC__) && (defined(__x86_64__) || defined(__i386__)) +/* GCC-compatible compiler, targeting x86/x86-64 */ +#include +#elif defined(__GNUC__) && defined(__ARM_NEON__) +/* GCC-compatible compiler, targeting ARM with NEON */ +#include +#elif defined(__GNUC__) && defined(__IWMMXT__) +/* GCC-compatible compiler, targeting ARM with WMMX */ + #include +#elif (defined(__GNUC__) || defined(__xlC__)) && (defined(__VEC__) || defined(__ALTIVEC__)) + /* XLC or GCC-compatible compiler, targeting PowerPC with VMX/VSX */ + #include +#elif defined(__GNUC__) && defined(__SPE__) + /* GCC-compatible compiler, targeting PowerPC with SPE */ + #include +#endif + +void update_progress(int32_t total_progress, int32_t page_progress, int32_t page_number, double eta) { + int extra_spaces = 2; + printf("\rProgress: %3d%% | Page %d progress: %3d%% | ETA: %.0fs%*s", total_progress, page_number, page_progress, eta, extra_spaces, ""); + printf("\033[%dD", extra_spaces); + fflush(stdout); +} + +void bgra_to_rgba(uint32_t *pixels, int tile_width, int tile_height) { + int num_pixels = tile_width * tile_height; + int num_pixels_aligned = (num_pixels / 4) * 4; + +#if defined(__ARM_NEON) + for (int i = 0; i < num_pixels_aligned; i += 4) { + uint32x4_t bgra = vld1q_u32(pixels + i); + uint32x4_t b_mask = vdupq_n_u32(0x000000FF); + uint32x4_t r_mask = vdupq_n_u32(0x00FF0000); + uint32x4_t b = vandq_u32(bgra, b_mask); + uint32x4_t r = vandq_u32(bgra, r_mask); + uint32x4_t br_swapped = vorrq_u32(vshlq_n_u32(b, 16), vshrq_n_u32(r, 16)); + uint32x4_t ga_alpha_mask = vdupq_n_u32(0xFF00FF00); + uint32x4_t ga_alpha = vandq_u32(bgra, ga_alpha_mask); + uint32x4_t rgba = vorrq_u32(ga_alpha, br_swapped); + vst1q_u32(pixels + i, rgba); + } +#elif defined(__SSE2__) + for (int i = 0; i < num_pixels_aligned; i += 4) { + __m128i bgra = _mm_loadu_si128((__m128i*)(pixels + i)); + __m128i b_mask = _mm_set1_epi32(0x000000FF); + __m128i r_mask = _mm_set1_epi32(0x00FF0000); + __m128i b = _mm_and_si128(bgra, b_mask); + __m128i r = _mm_and_si128(bgra, r_mask); + __m128i br_swapped = _mm_or_si128(_mm_slli_epi32(b, 16), _mm_srli_epi32(r, 16)); + __m128i ga_alpha_mask = _mm_set1_epi32(0xFF00FF00); + __m128i ga_alpha = _mm_and_si128(bgra, ga_alpha_mask); + __m128i rgba = _mm_or_si128(ga_alpha, br_swapped); + _mm_storeu_si128((__m128i*)(pixels + i), rgba); + } +#else + for (int i = num_pixels_aligned; i < num_pixels; ++i) { + uint32_t val = pixels[i]; + pixels[i] = ((val & 0xff) << 16) | (val & 0x00ff00) | ((val & 0xff0000) >> 16) | (val & 0xff000000); + } +#endif +} + +void write_page_to_tiff(TIFF *output_tiff, isyntax_t *isyntax, isyntax_cache_t *isyntax_cache, isyntax_level_t *level, + int32_t tile_width, int32_t tile_height, int32_t *total_tiles_written, int32_t total_tiles, + clock_t global_start_time, uint16_t compression_type, uint16_t quality) { + int32_t width = libisyntax_level_get_width(level); + int32_t height = libisyntax_level_get_height(level); + int32_t scale = libisyntax_level_get_scale(level); + + // Set the TIFF properties for the current level. + TIFFSetField(output_tiff, TIFFTAG_IMAGEWIDTH, width); + TIFFSetField(output_tiff, TIFFTAG_IMAGELENGTH, height); + TIFFSetField(output_tiff, TIFFTAG_BITSPERSAMPLE, 8); + TIFFSetField(output_tiff, TIFFTAG_SAMPLESPERPIXEL, 4); + + if (compression_type == COMPRESSION_JPEG) { + TIFFSetField(output_tiff, TIFFTAG_COMPRESSION, COMPRESSION_JPEG); + TIFFSetField(output_tiff, TIFFTAG_JPEGQUALITY, quality); + } else if (compression_type == COMPRESSION_LZW) { + TIFFSetField(output_tiff, TIFFTAG_COMPRESSION, COMPRESSION_LZW); + } + + TIFFSetField(output_tiff, TIFFTAG_PHOTOMETRIC, PHOTOMETRIC_RGB); + TIFFSetField(output_tiff, TIFFTAG_ORIENTATION, ORIENTATION_TOPLEFT); + TIFFSetField(output_tiff, TIFFTAG_PLANARCONFIG, PLANARCONFIG_CONTIG); + TIFFSetField(output_tiff, TIFFTAG_TILEWIDTH, tile_width); + TIFFSetField(output_tiff, TIFFTAG_TILELENGTH, tile_height); + + // Set the resolution + double res_x = 10000.0 / libisyntax_level_get_mpp_x(level); + TIFFSetField(output_tiff, TIFFTAG_XRESOLUTION, res_x); + double res_y = 10000.0 / libisyntax_level_get_mpp_y(level); + TIFFSetField(output_tiff, TIFFTAG_YRESOLUTION, res_y); + TIFFSetField(output_tiff, TIFFTAG_RESOLUTIONUNIT, RESUNIT_CENTIMETER); + TIFFSetField(output_tiff, TIFFTAG_EXTRASAMPLES, 1, (uint16_t[]) {EXTRASAMPLE_UNASSALPHA}); + + if (scale == 0) { + TIFFSetField(output_tiff, TIFFTAG_SUBFILETYPE, 0); + } else { + TIFFSetField(output_tiff, TIFFTAG_SUBFILETYPE, FILETYPE_REDUCEDIMAGE); + } + + int32_t tile_progress = 0; + int32_t tiles_in_page = ((height + tile_height - 1) / tile_height) * ((width + tile_width - 1) / tile_width); + + for (int32_t y_coord = 0; y_coord < height; y_coord += tile_height) { + for (int32_t x_coord = 0; x_coord < width; x_coord += tile_width) { + // At the borders the tile size can be smaller + int32_t region_width = (x_coord + tile_width > width) ? width - x_coord : tile_width; + int32_t region_height = (y_coord + tile_height > height) ? height - y_coord : tile_height; + + uint32_t *pixels = NULL; + + // Checking if we are reading within bounds + assert(x_coord + region_width <= width); + assert(y_coord + region_height <= height); + + assert(libisyntax_read_region(isyntax, isyntax_cache, scale, x_coord, y_coord, region_width, region_height, + &pixels) == LIBISYNTAX_OK); + + // Convert data to the correct pixel format (bgra->rgba). + bgra_to_rgba(pixels, region_width, region_height); + + uint32_t *tile_pixels = NULL; + + // In case our actual tile is smaller, we need to convert it to a full tile. + if (region_width != tile_width || region_height != tile_height) { + tile_pixels = calloc(tile_width * tile_height, sizeof(uint32_t)); + for (int32_t row = 0; row < region_height; ++row) { + memcpy(tile_pixels + row * tile_width, pixels + row * region_width, + region_width * sizeof(uint32_t)); + } + libisyntax_tile_free_pixels(pixels); // Move this line here + } else { + tile_pixels = pixels; + } + + // Write the tile to the output TIFF. + TIFFWriteTile(output_tiff, tile_pixels, x_coord, y_coord, 0, 0); + if (tile_pixels != NULL) { + libisyntax_tile_free_pixels(tile_pixels); + } + + ++tile_progress; + *total_tiles_written += 1; + int32_t tile_percent = (tile_progress * 100) / tiles_in_page; + int32_t total_progress = ((*total_tiles_written) * 100) / total_tiles; + + clock_t current_global_time = clock(); + double elapsed_global_time = (double) (current_global_time - global_start_time) / CLOCKS_PER_SEC; + double avg_time_per_tile = (*total_tiles_written > 0) ? (elapsed_global_time / *total_tiles_written) : 0; + double eta = avg_time_per_tile * (total_tiles - *total_tiles_written); + update_progress(total_progress, tile_percent, scale, eta); + } + } + + // Write the directory for the current level. + TIFFWriteDirectory(output_tiff); +} + +uint64_t parse_cache_size(const char *size_str) { + uint64_t size; + char unit; + + if (sscanf(size_str, "%lld%c", &size, &unit) == 2) { + if (unit == 'M') { + if (size > INT64_MAX / (1024)) { + printf("Error: Cache size too large.\n"); + return -1; + } + size *= 1024; + } else if (unit == 'G') { + if (size > INT64_MAX / (1024 * 1024)) { + printf("Error: Cache size too large.\n"); + return -1; + } + size *= 1024 * 1024; + } else { + printf("Error: Invalid unit for cache size. Use 'M' for megabytes or 'G' for gigabytes.\n"); + return -1; + } + } else if (sscanf(size_str, "%lld", &size) != 1) { + printf("Error: Invalid cache size format.\n"); + return -1; + } + + return size; +} + + +int main(int argc, char **argv) { + const char *usage_string = + "Usage: isyntax-to-tiff [OPTIONS] INPUT OUTPUT\n\n" + "Converts Philips iSyntax files to a multi-resolution BigTIFF file.\n\n" + "Positional arguments:\n" + " INPUT Path to the input iSyntax file.\n" + " OUTPUT Path to the output TIFF file.\n\n" + "Options:\n" + " --help Show this help message and exit.\n\n" + " --start-at-page PAGE Specifies the page to start at (default: 0).\n" + " --tile-size SIZE Specifies the tile size for the output TIFF (default: 1024).\n" + " Must be a positive integer.\n\n" + " --compression TYPE Specifies the compression type for the output TIFF.\n" + " Supported types: JPEG, LZW, NONE (default: JPEG).\n\n" + " --quality VALUE Specifies the quality for JPEG compression (0-100).\n" + " Only applicable when using JPEG compression (default: 80).\n\n" + " --cache-size SIZE Specifies the cache size for the iSyntax library.\n" + " Accepts a number followed by 'M' (for megabytes) or 'G' (for gigabytes),\n" + " or just a number for kilobytes (default: 2000).\n\n" + "Example:\n\n" + " isyntax-to-tiff --tile-size 512 --compression JPEG --quality 90 --cache-size 1G input.isyntax output.tiff\n\n" + "This command will convert the input.isyntax file into an output.tiff file with a tile size of 512, JPEG compression at 90 quality, and a cache size of 1 gigabyte.\n"; + + if (argc < 3) { + printf("Error: Missing input and/or output file arguments.\n\n"); + printf("%s", usage_string); + return -1; + } + + char *filename = argv[1]; + char *output_tiffname = argv[2]; + + uint64_t cache_size = 2000; + int32_t tile_size = 1024; + int start_at_page = 0; + + int compression_type = COMPRESSION_JPEG; + int quality = 80; + + for (int i = 3; i < argc; ++i) { + if (strcmp(argv[i], "--tile-size") == 0) { + if (i + 1 < argc) { + tile_size = atoi(argv[i + 1]); + if (tile_size <= 0) { + printf("Error: Invalid tile size. Please provide a positive integer value for the tile size.\n"); + return -1; + } + i++; // Skip the next argument (tile size value) + } else { + printf("Error: Missing value for --tile-size option.\n"); + return -1; + } + } else if (strcmp(argv[i], "--compression") == 0) { + if (i + 1 < argc) { + if (strcmp(argv[i + 1], "JPEG") == 0) { + compression_type = COMPRESSION_JPEG; + } else if (strcmp(argv[i + 1], "LZW") == 0) { + compression_type = COMPRESSION_LZW; + } else if (strcmp(argv[i + 1], "NONE") == 0) { + compression_type = COMPRESSION_NONE; + } else { + printf("Error: Invalid compression type. Supported types are JPEG, LZW, and NONE.\n"); + return -1; + } + i++; // Skip the next argument (compression type value) + } else { + printf("Error: Missing value for --compression option.\n"); + return -1; + } + } else if (strcmp(argv[i], "--quality") == 0) { + if (i + 1 < argc) { + quality = atoi(argv[i + 1]); + if (quality < 0 || quality > 100) { + printf("Error: Invalid quality value. Please provide an integer value between 0 and 100 for the quality.\n"); + return -1; + } + if (compression_type != COMPRESSION_JPEG) { + printf("Warning: The --quality flag is ignored with the current compression type. Quality is only applicable to JPEG compressions.\n"); + } + i++; // Skip the next argument (quality value) + } else { + printf("Error: Missing value for --quality option.\n"); + return -1; + } + } else if (strcmp(argv[i], "--help") == 0) { + printf("%s", usage_string); + return 0; + } else if (strcmp(argv[i], "--start-at-page") == 0) { + if (i + 1 < argc) { + start_at_page = atoi(argv[i + 1]); + if (start_at_page < 0) { + printf("Error: Invalid page number. Please provide a positive integer value for the page number.\n"); + return -1; + } + i++; // Skip the next argument (page number value) + } else { + printf("Error: Missing value for --start-at-page option.\n"); + return -1; + } + } else if (strcmp(argv[i], "--cache-size") == 0) { + if (i + 1 < argc) { + cache_size = parse_cache_size(argv[i + 1]); + if (cache_size >= INT64_MAX || cache_size < 0) { + printf("Error: Cache size not suitable for the system.\n"); + return -1; + } + i++; // Skip the next argument (cache size value) + } else { + printf("Error: Missing value for --cache-size option.\n"); + return -1; + } + + } else { + printf("Error: Unknown option %s\n", argv[i]); + return -1; + } + } + int32_t tile_width = tile_size; + int32_t tile_height = tile_size; + + libisyntax_init(); + + isyntax_t *isyntax; + if (libisyntax_open(filename, /*is_init_allocators=*/0, &isyntax) != LIBISYNTAX_OK) { + fprintf(stderr, "Failed to open %s\n", filename); + return -1; + } + + + isyntax_cache_t *isyntax_cache = NULL; + if (libisyntax_cache_create("isyntax-to-tiff cache", cache_size, &isyntax_cache) != LIBISYNTAX_OK) { + fprintf(stderr, "Failed to create iSyntax cache with size %llu.\n", cache_size); + libisyntax_close(isyntax); + return -1; + } + if (libisyntax_cache_inject(isyntax_cache, isyntax) != LIBISYNTAX_OK) { + fprintf(stderr, "Failed to inject iSyntax cache into iSyntax instance.\n"); + libisyntax_cache_destroy(isyntax_cache); + libisyntax_close(isyntax); + return -1; + } + + // Initialize the output TIFF file. + TIFF *output_tiff; + output_tiff = TIFFOpen(output_tiffname, "w8"); + if (!output_tiff) { + fprintf(stderr, "Failed to create %s\n", output_tiffname); + return -1; + } + + const isyntax_image_t *image = libisyntax_get_image(isyntax, 0); + + // Check if the page we selected actually exists. + int32_t level_count = libisyntax_image_get_level_count(image); + if (start_at_page >= level_count) { + fprintf(stderr, "Error: The page number %d is out of range. The image only has %d pages. Set --start-at-page to a smaller value.\n", start_at_page, level_count); + libisyntax_cache_destroy(isyntax_cache); + libisyntax_close(isyntax); + TIFFClose(output_tiff); + return -1; + } + + int32_t num_levels = libisyntax_image_get_level_count(image); + int32_t total_tiles = 0; + + // Let's find the total number of tiles so we can have a progress counter + for (int32_t level = start_at_page; level < num_levels; ++level) { + isyntax_level_t *current_level = libisyntax_image_get_level(image, level); + int32_t width = libisyntax_level_get_width(current_level); + int32_t height = libisyntax_level_get_height(current_level); + int32_t tiles_in_page = ((height + tile_height - 1) / tile_height) * ((width + tile_width - 1) / tile_width); + + total_tiles += tiles_in_page; + } + + int32_t total_tiles_written = 0; + clock_t global_start_time = clock(); + for (int32_t level = start_at_page; level < num_levels; ++level) { + isyntax_level_t *current_level = libisyntax_image_get_level(image, level); + write_page_to_tiff(output_tiff, isyntax, isyntax_cache, current_level, tile_width, tile_height, + &total_tiles_written, total_tiles, global_start_time, compression_type, quality); + } + + // Close the output TIFF file. + TIFFClose(output_tiff); + + // Clean up. + libisyntax_cache_destroy(isyntax_cache); + libisyntax_close(isyntax); + + return 0; +} diff --git a/src/libisyntax.c b/src/libisyntax.c index eb90f70..00b6c99 100644 --- a/src/libisyntax.c +++ b/src/libisyntax.c @@ -28,6 +28,7 @@ #define STB_SPRINTF_IMPLEMENTATION #include "stb_sprintf.h" + #define STB_IMAGE_IMPLEMENTATION #include "stb_image.h" @@ -38,6 +39,7 @@ #include "libisyntax.h" #include "isyntax.h" #include "isyntax_reader.h" +#include static platform_thread_info_t thread_infos[MAX_THREAD_COUNT]; @@ -310,6 +312,94 @@ isyntax_error_t libisyntax_tile_read(isyntax_t* isyntax, isyntax_cache_t* isynta return LIBISYNTAX_OK; } +isyntax_error_t libisyntax_read_region(isyntax_t* isyntax, isyntax_cache_t* isyntax_cache, int32_t level, + int64_t x, int64_t y, int64_t width, int64_t height, uint32_t** out_pixels) { + + // Get the level + assert(level < isyntax->images[0].level_count); + isyntax_level_t* current_level = &isyntax->images[0].levels[level]; + + int32_t PER_LEVEL_PADDING = 3; + int32_t num_levels = isyntax->images[0].level_count; + int32_t offset = ((PER_LEVEL_PADDING << num_levels) - PER_LEVEL_PADDING) >> level; + + x += offset; + y += offset; + + int32_t tile_width = isyntax->tile_width; + int32_t tile_height = isyntax->tile_height; + + int64_t start_tile_x = x / tile_width; + int64_t end_tile_x = (x + width - 1) / tile_width; + int64_t start_tile_y = y / tile_height; + int64_t end_tile_y = (y + height - 1) / tile_height; + + // Allocate memory for region + *out_pixels = (uint32_t*)malloc(width * height * sizeof(uint32_t)); + + // Initialize the empty tile as a NULL pointer + // TODO: Maybe you want to attach this to some object so we do not need to keep on reallocating these things? + uint32_t *empty_tile = NULL; + + // Read tiles and copy the relevant portion of each tile to the region + for (int64_t tile_y = start_tile_y; tile_y <= end_tile_y; ++tile_y) { + for (int64_t tile_x = start_tile_x; tile_x <= end_tile_x; ++tile_x) { + // Calculate the portion of the tile to be copied + int64_t src_x = (tile_x == start_tile_x) ? x % tile_width : 0; + int64_t src_y = (tile_y == start_tile_y) ? y % tile_height : 0; + int64_t dest_x = (tile_x == start_tile_x) ? 0 : (tile_x - start_tile_x) * tile_width - (x % tile_width); + int64_t dest_y = (tile_y == start_tile_y) ? 0 : (tile_y - start_tile_y) * tile_height - (y % tile_height); + int64_t copy_width = (tile_x == end_tile_x) ? (x + width - 1) % tile_width - src_x + 1 : tile_width - src_x; + int64_t copy_height = (tile_y == end_tile_y) ? (y + height - 1) % tile_height - src_y + 1 : tile_height - src_y; + + + uint32_t *pixels = NULL; + + int64_t tile_index = tile_y * current_level->width_in_tiles + tile_x; + // Check if tile exists, if not, don't use the function to read the tile and immediately return an empty + // tile. + bool tile_exists = (isyntax->images[0].levels[level].tiles + tile_index)->exists; + if (tile_exists) { + // Read tile + assert(tile_x < current_level->width_in_tiles); + assert(tile_y < current_level->height_in_tiles); + assert(libisyntax_tile_read(isyntax, isyntax_cache, level, tile_x, tile_y, &pixels) == LIBISYNTAX_OK); + } else { + // Allocate memory for the empty tile and fill it with non-transparent white pixels only when required + if (empty_tile == NULL) { + empty_tile = (uint32_t*)malloc(tile_width * tile_height * sizeof(uint32_t)); + for (int64_t i = 0; i < tile_height; ++i) { + for (int64_t j = 0; j < tile_width; ++j) { + empty_tile[i * tile_width + j] = 0xFFFFFFFFu; // Could be 0x00FFFFFFu for A=0 + } + } + } + pixels = empty_tile; + } + + // Copy the relevant portion of the tile to the region + for (int64_t i = 0; i < copy_height; ++i) { + int64_t dest_index = (dest_y + i) * width + dest_x; + int64_t src_index = (src_y + i) * tile_width + src_x; + memcpy((*out_pixels) + dest_index, + pixels + src_index, + copy_width * sizeof(uint32_t)); + } + // Free the tile data if it exists + if (pixels != NULL && pixels != empty_tile) { + libisyntax_tile_free_pixels(pixels); + } + } + } + + // Free the empty tile data if it was allocated + if (empty_tile != NULL) { + libisyntax_tile_free_pixels(empty_tile); + } + + return LIBISYNTAX_OK; +} + void libisyntax_tile_free_pixels(uint32_t* pixels) { free(pixels); } diff --git a/src/libisyntax.h b/src/libisyntax.h index 4fbdb0b..ab24618 100644 --- a/src/libisyntax.h +++ b/src/libisyntax.h @@ -70,6 +70,8 @@ void libisyntax_cache_destroy(isyntax_cache_t* isyntax_cache); // allow user to supply the buffer, as long as the user can compute the size. isyntax_error_t libisyntax_tile_read(isyntax_t* isyntax, isyntax_cache_t* isyntax_cache, int32_t level, int64_t tile_x, int64_t tile_y, uint32_t** out_pixels); +isyntax_error_t libisyntax_read_region(isyntax_t* isyntax, isyntax_cache_t* isyntax_cache, int32_t level, + int64_t x, int64_t y, int64_t width, int64_t height, uint32_t** out_pixels); void libisyntax_tile_free_pixels(uint32_t* pixels);