commit acb5335cb8b815da1518b97f8b1cdd23157f4914
parent d5a74b886ef518818298582d821ff2aaceff8df8
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Tue, 15 May 2018 11:15:59 +0200
Rename the htcop project in htcp
Improve the README and refactor the CMakeLists files
Diffstat:
23 files changed, 1774 insertions(+), 1718 deletions(-)
diff --git a/README.md b/README.md
@@ -0,0 +1,35 @@
+# High-Tune: Cloud Properties
+
+This project describes the htcp binary fileformat that is used to store the
+some properties of a cloud. The provided `les2htcp` command line tool extracts
+these properties from a LES output stored in a
+[NetCDF](https://www.unidata.ucar.edu/software/netcdf/) file, and converts them
+in the `htcp` fileformat. Finally, a library is provided to load the `htcp`
+fileformat. Since the data to load can be huge, this library silently
+loads/unloads them dynamically allowing to process data that are too large to
+fit in the main memory.
+
+## How to build
+
+The `les2htcp` program and `htcp` library are compatible GNU/Linux 64-bits.
+They rely on the [CMake](http://www.cmake.org) and the
+[RCMake](https://gitlab.com/vaplv/rcmake/) packages to build. They also depend
+on the [RSys](https://gitlab.com/vaplv/rsys/) library. Furthermore, the
+`les2htcp` tool depends on the
+[NetCDF](https://www.unidata.ucar.edu/software/netcdf/) C library.
+
+To build them, first ensure that CMake is installed on your system. Then
+install the RCMake package as well as the aforementioned prerequisites. Finally
+generate the project from the `cmake/CMakeLists.txt` file by appending to the
+`CMAKE_PREFIX_PATH` variable the install directories of its dependencies. The
+resulting project can be edited, built, tested and installed as any CMake
+project. Refer to the [CMake](https://cmake.org/documentation) for further
+informations on CMake.
+
+## Licenses
+
+The `les2htcp` tool and the `htcp` library are free software copyright
+(C) |Meso|Star 2018. They are released under the GPL v3+ license: GNU GPL
+version 3 or later. You are welcome to redistribute them under certain
+conditions; refer to the COPYING file for details.
+
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -14,21 +14,14 @@
# along with this program. If not, see <http://www.gnu.org/licenses/>.
cmake_minimum_required(VERSION 2.8)
-project(htcop C)
+project(htcp)
enable_testing()
option(NO_TEST "Do not build tests" OFF)
-add_subdirectory(les2htcop) # Build utility
-
-set(HTCOP_SOURCE_DIR ${PROJECT_SOURCE_DIR}/../src)
-
################################################################################
-# Check dependencies
+# Prerequisites
################################################################################
-get_filename_component(_current_source_dir ${CMAKE_CURRENT_LIST_FILE} PATH)
-set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${_current_source_dir}/)
-
find_package(RCMake 0.3 REQUIRED)
find_package(RSys 0.6 REQUIRED)
@@ -39,62 +32,24 @@ include(rcmake_runtime)
include_directories(${RSys_INCLUDE_DIR})
################################################################################
-# Configure and define targets
-################################################################################
-set(VERSION_MAJOR 0)
-set(VERSION_MINOR 0)
-set(VERSION_PATCH 0)
-set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
-
-set(HTCOP_FILES_SRC htcop.c)
-set(HTCOP_FILES_INC )
-set(HTCOP_FILES_INC_API htcop.h)
-set(HTCOP_FILES_DOC COPYING)
-
-# Prepend each file in the `HTCOP_FILES_<SRC|INC>' list by `HTCOP_SOURCE_DIR'
-rcmake_prepend_path(HTCOP_FILES_SRC ${HTCOP_SOURCE_DIR})
-rcmake_prepend_path(HTCOP_FILES_INC ${HTCOP_SOURCE_DIR})
-rcmake_prepend_path(HTCOP_FILES_INC_API ${HTCOP_SOURCE_DIR})
-rcmake_prepend_path(HTCOP_FILES_DOC ${PROJECT_SOURCE_DIR}/../)
-
-add_library(htcop SHARED ${HTCOP_FILES_SRC} ${HTCOP_FILES_INC} ${HTCOP_FILES_INC_API})
-target_link_libraries(htcop RSys)
-
-set_target_properties(htcop PROPERTIES
- DEFINE_SYMBOL HTCOP_SHARED_BUILD
- VERSION ${VERSION}
- SOVERSION ${VERSION_MAJOR})
-
-rcmake_setup_devel(htcop HTCOP ${VERSION} high_tune/htcop.h)
-
-################################################################################
-# Add tests
+# Test utilities
################################################################################
if(NOT NO_TEST)
function(build_test _name)
add_executable(${_name}
- ${HTCOP_SOURCE_DIR}/${_name}.c)
- target_link_libraries(${_name} htcop)
+ ${HTCP_SOURCE_DIR}/${_name}.c)
+ target_link_libraries(${_name} htcp)
endfunction()
function(new_test _name)
build_test(${_name})
add_test(${_name} ${_name})
endfunction()
-
- new_test(test_htcop)
- new_test(test_htcop_load)
-
- build_test(test_htcop_load_from_file)
endif()
################################################################################
-# Define output & install directories
+# Sub projects
################################################################################
-install(TARGETS htcop
- ARCHIVE DESTINATION bin
- LIBRARY DESTINATION lib
- RUNTIME DESTINATION bin)
-install(FILES ${HTCOP_FILES_INC_API} DESTINATION include/high_tune)
-install(FILES ${HTCOP_FILES_DOC} DESTINATION share/doc/htcop)
+add_subdirectory(htcp)
+add_subdirectory(les2htcp)
diff --git a/cmake/htcp/CMakeLists.txt b/cmake/htcp/CMakeLists.txt
@@ -0,0 +1,68 @@
+# Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+cmake_minimum_required(VERSION 2.8)
+project(htcp-library C)
+
+set(HTCP_SOURCE_DIR ${PROJECT_SOURCE_DIR}/../../src)
+
+################################################################################
+# Configure and define targets
+################################################################################
+set(VERSION_MAJOR 0)
+set(VERSION_MINOR 0)
+set(VERSION_PATCH 0)
+set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
+
+set(HTCP_FILES_SRC htcp.c)
+set(HTCP_FILES_INC )
+set(HTCP_FILES_INC_API htcp.h)
+set(HTCP_FILES_DOC COPYING)
+
+# Prepend each file in the `HTCP_FILES_<SRC|INC>' list by `HTCP_SOURCE_DIR'
+rcmake_prepend_path(HTCP_FILES_SRC ${HTCP_SOURCE_DIR})
+rcmake_prepend_path(HTCP_FILES_INC ${HTCP_SOURCE_DIR})
+rcmake_prepend_path(HTCP_FILES_INC_API ${HTCP_SOURCE_DIR})
+rcmake_prepend_path(HTCP_FILES_DOC ${PROJECT_SOURCE_DIR}/../)
+
+add_library(htcp SHARED ${HTCP_FILES_SRC} ${HTCP_FILES_INC} ${HTCP_FILES_INC_API})
+target_link_libraries(htcp RSys)
+
+set_target_properties(htcp PROPERTIES
+ DEFINE_SYMBOL HTCP_SHARED_BUILD
+ VERSION ${VERSION}
+ SOVERSION ${VERSION_MAJOR})
+
+rcmake_setup_devel(htcp HTCP ${VERSION} high_tune/htcp.h)
+
+################################################################################
+# Add tests
+################################################################################
+if(NOT NO_TEST)
+ new_test(test_htcp)
+ new_test(test_htcp_load)
+ build_test(test_htcp_load_from_file)
+endif()
+
+################################################################################
+# Define output & install directories
+################################################################################
+install(TARGETS htcp
+ ARCHIVE DESTINATION bin
+ LIBRARY DESTINATION lib
+ RUNTIME DESTINATION bin)
+install(FILES ${HTCP_FILES_INC_API} DESTINATION include/high_tune)
+install(FILES ${HTCP_FILES_DOC} DESTINATION share/doc/htcp)
+
diff --git a/cmake/les2htcop/CMakeLists.txt b/cmake/les2htcop/CMakeLists.txt
@@ -1,78 +0,0 @@
-# Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with this program. If not, see <http://www.gnu.org/licenses/>.
-
-cmake_minimum_required(VERSION 2.8)
-project(les2htcop C)
-enable_testing()
-
-set(LES2HTCOP_SOURCE_DIR ${PROJECT_SOURCE_DIR}/../../src)
-option(NO_TEST "Do not build tests" OFF)
-
-################################################################################
-# Check dependencies
-################################################################################
-get_filename_component(_current_source_dir ${CMAKE_CURRENT_LIST_FILE} PATH)
-set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${_current_source_dir}/../)
-
-find_package(RCMake 0.3 REQUIRED)
-find_package(RSys 0.6 REQUIRED)
-find_package(NetCDF REQUIRED)
-
-set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${RCMAKE_SOURCE_DIR})
-include(rcmake)
-include(rcmake_runtime)
-
-include_directories(
- ${RSys_INCLUDE_DIR}
- ${NETCDF_C_INCLUDE_DIRS}
- ${CMAKE_CURRENT_BINARY_DIR})
-
-################################################################################
-# Generate files
-################################################################################
-set(VERSION_MAJOR 0)
-set(VERSION_MINOR 0)
-set(VERSION_PATCH 0)
-set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
-
-configure_file(${LES2HTCOP_SOURCE_DIR}/les2htcop.h.in
- ${CMAKE_CURRENT_BINARY_DIR}/les2htcop.h @ONLY)
-
-################################################################################
-# Configure and define targets
-################################################################################
-set(LES2HTCOP_FILES_SRC les2htcop.c)
-set(LES2HTCOP_FILES_INC les2htcop.h.in)
-set(LES2HTCOP_FILES_DOC COPYING)
-
-# Prepend each file in the `HTNC_FILES_<SRC|INC>' list by `LES2HTCOP_SOURCE_DIR'
-rcmake_prepend_path(LES2HTCOP_FILES_SRC ${LES2HTCOP_SOURCE_DIR})
-rcmake_prepend_path(LES2HTCOP_FILES_INC ${LES2HTCOP_SOURCE_DIR})
-rcmake_prepend_path(LES2HTCOP_FILES_DOC ${PROJECT_SOURCE_DIR}/../../)
-
-add_executable(les2htcop ${LES2HTCOP_FILES_SRC})
-target_link_libraries(les2htcop RSys ${NETCDF_C_LIBRARIES})
-
-set_target_properties(les2htcop PROPERTIES
- DEFINE_SYMBOL LES2HTCOP_SHARED_BUILD
- VERSION ${VERSION}
- SOVERSION ${VERSION_MAJOR})
-
-################################################################################
-# Define output & install directories
-################################################################################
-install(TARGETS les2htcop RUNTIME DESTINATION bin)
-install(FILES ${LES2HTCOP_FILES_DOC} DESTINATION share/doc/les2htcop)
-
diff --git a/cmake/les2htcp/CMakeLists.txt b/cmake/les2htcp/CMakeLists.txt
@@ -0,0 +1,69 @@
+# Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+cmake_minimum_required(VERSION 2.8)
+project(les2htcp C)
+
+set(LES2HTCP_SOURCE_DIR ${PROJECT_SOURCE_DIR}/../../src)
+
+################################################################################
+# Check dependencies
+################################################################################
+get_filename_component(_current_source_dir ${CMAKE_CURRENT_LIST_FILE} PATH)
+set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${_current_source_dir})
+
+find_package(NetCDF REQUIRED)
+
+include_directories(
+ ${NETCDF_C_INCLUDE_DIRS}
+ ${CMAKE_CURRENT_BINARY_DIR})
+
+################################################################################
+# Generate files
+################################################################################
+set(VERSION_MAJOR 0)
+set(VERSION_MINOR 0)
+set(VERSION_PATCH 0)
+set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
+
+configure_file(${LES2HTCP_SOURCE_DIR}/les2htcp.h.in
+ ${CMAKE_CURRENT_BINARY_DIR}/les2htcp.h @ONLY)
+
+################################################################################
+# Configure and define targets
+################################################################################
+set(LES2HTCP_FILES_SRC les2htcp.c)
+set(LES2HTCP_FILES_INC les2htcp.h.in)
+set(LES2HTCP_FILES_DOC COPYING)
+
+# Prepend each file in the `HTNC_FILES_<SRC|INC>' list by `LES2HTCP_SOURCE_DIR'
+rcmake_prepend_path(LES2HTCP_FILES_SRC ${LES2HTCP_SOURCE_DIR})
+rcmake_prepend_path(LES2HTCP_FILES_INC ${LES2HTCP_SOURCE_DIR})
+rcmake_prepend_path(LES2HTCP_FILES_DOC ${PROJECT_SOURCE_DIR}/../../)
+
+add_executable(les2htcp ${LES2HTCP_FILES_SRC})
+target_link_libraries(les2htcp RSys ${NETCDF_C_LIBRARIES})
+
+set_target_properties(les2htcp PROPERTIES
+ DEFINE_SYMBOL LES2HTCP_SHARED_BUILD
+ VERSION ${VERSION}
+ SOVERSION ${VERSION_MAJOR})
+
+################################################################################
+# Define output & install directories
+################################################################################
+install(TARGETS les2htcp RUNTIME DESTINATION bin)
+install(FILES ${LES2HTCP_FILES_DOC} DESTINATION share/doc/les2htcp)
+
diff --git a/cmake/FindNetCDF.cmake b/cmake/les2htcp/FindNetCDF.cmake
diff --git a/doc/htcop.txt b/doc/htcop.txt
@@ -1,29 +0,0 @@
-<htnc> ::= <header>
- <definition>
- <lower-pos>
- <voxel-size>
- <padding>
- <RVT> # align on pagesize
- <padding>
- <RCT> # align on pagesize
- <padding>
-
-<header> ::= <pagesize> <is-Z-irregular>
-<definition> ::= <#X> <#Y> <#Z> <#time>
-<lower-pos> ::= <double3>
-<voxel-size> ::= <double3> [DOUBLE ... ] # optionnal double for irregular Z
-<RVT> ::= DOUBLE [DOUBLE ... ]
-<RCT> ::= DOUBLE [DOUBLE ... ]
-
-<double3> ::= DOUBLE DOUBLE DOUBLE
-
-<#X> ::= INT32
-<#Y> ::= INT32
-<#Z> ::= INT32
-<#time> ::= INT32
-
-<is-Z-irregular> ::= INT8
-
-<pagesize> ::= INT64
-
-<padding> ::= [ BYTES ... ] # Alignment address onto <pagesize>
diff --git a/doc/htcp.txt b/doc/htcp.txt
@@ -0,0 +1,36 @@
+# Copyright (C) 2018 |Meso|Star> <contact@meso-star.com>
+#
+# Copying and distribution of this file, with or without modification,
+# are permitted in any medium without royalty provided the copyright
+# notice and this notice are preserved. This file is offered as-is,
+# without any warranty.
+
+<htnc> ::= <header>
+ <definition>
+ <lower-pos>
+ <voxel-size>
+ <padding>
+ <RVT> # align on pagesize
+ <padding>
+ <RCT> # align on pagesize
+ <padding>
+
+<header> ::= <pagesize> <is-Z-irregular>
+<definition> ::= <#X> <#Y> <#Z> <#time>
+<lower-pos> ::= <double3>
+<voxel-size> ::= <double3> [DOUBLE ... ] # optionnal double for irregular Z
+<RVT> ::= DOUBLE [DOUBLE ... ]
+<RCT> ::= DOUBLE [DOUBLE ... ]
+
+<double3> ::= DOUBLE DOUBLE DOUBLE
+
+<#X> ::= INT32
+<#Y> ::= INT32
+<#Z> ::= INT32
+<#time> ::= INT32
+
+<is-Z-irregular> ::= INT8
+
+<pagesize> ::= INT64
+
+<padding> ::= [ BYTES ... ] # Alignment address onto <pagesize>
diff --git a/src/htcop.c b/src/htcop.c
@@ -1,335 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#define _POSIX_C_SOURCE 200809L /* mmap support */
-#define _DEFAULT_SOURCE 1 /* MAP_POPULATE support */
-
-#include "htcop.h"
-
-#include <rsys/dynamic_array_double.h>
-#include <rsys/logger.h>
-#include <rsys/ref_count.h>
-#include <rsys/mem_allocator.h>
-
-#include <errno.h>
-#include <string.h>
-#include <sys/mman.h>
-#include <unistd.h>
-
-enum { X, Y, Z, TIME }; /* Helper constants */
-
-struct htcop {
- int64_t pagesize;
- int8_t irregular_z;
- int32_t definition[4];
- double lower[3];
- double vxsz[2]; /* Size of the voxels in X and Y */
- struct darray_double vxsz_z; /* Size of the voxels along the Z dimension */
-
- double* RCT; /* Mapped memory */
- double* RVT; /* Mapped memory */
- size_t RCT_length; /* In bytes */
- size_t RVT_length; /* In bytes */
-
- size_t pagesize_os; /* Page size of the os */
- int verbose; /* Verbosity level */
- struct logger* logger;
- struct mem_allocator* allocator;
- ref_T ref;
-};
-
-/*******************************************************************************
- * Local functions
- ******************************************************************************/
-static void
-log_msg
- (const struct htcop* htcop,
- const enum log_type stream,
- const char* msg,
- va_list vargs)
-{
- ASSERT(htcop && msg);
- if(htcop->verbose) {
- res_T res; (void)res;
- res = logger_vprint(htcop->logger, stream, msg, vargs);
- ASSERT(res == RES_OK);
- }
-}
-
-static void
-log_err(const struct htcop* htcop, const char* msg, ...)
-{
- va_list vargs_list;
- ASSERT(htcop && msg);
- va_start(vargs_list, msg);
- log_msg(htcop, LOG_ERROR, msg, vargs_list);
- va_end(vargs_list);
-}
-
-static void
-reset_htcop(struct htcop* htcop)
-{
- htcop->pagesize = 0;
- htcop->irregular_z = 0;
- htcop->definition[0] = 0;
- htcop->definition[1] = 0;
- htcop->definition[2] = 0;
- htcop->lower[0] = -1;
- htcop->lower[1] = -1;
- htcop->lower[2] = -1;
- htcop->vxsz[0] = -1;
- htcop->vxsz[1] = -1;
- darray_double_clear(&htcop->vxsz_z);
- if(htcop->RCT) munmap(htcop->RCT, htcop->RCT_length);
- if(htcop->RVT) munmap(htcop->RVT, htcop->RVT_length);
- htcop->RCT = NULL;
- htcop->RVT = NULL;
- htcop->RCT_length = 0;
- htcop->RVT_length = 0;
-}
-
-static res_T
-load_stream(struct htcop* htcop, FILE* stream, const char* stream_name)
-{
- size_t nz = 0;
- size_t map_len = 0;
- off_t offset = 0;
- res_T res = RES_OK;
- ASSERT(htcop && stream && stream_name);
-
- reset_htcop(htcop);
-
- #define READ(Var, N, Name) { \
- if(fread((Var), sizeof(*(Var)), (N), stream) != (N)) { \
- fprintf(stderr, "%s: could not read the %s\n", stream_name, Name); \
- res = RES_IO_ERR; \
- goto error; \
- } \
- } (void)0
- READ(&htcop->pagesize, 1, "page size");
- if(!IS_POW2(htcop->pagesize)) {
- log_err(htcop, "%s: invalid page size `%li'. It must be a power of two.\n",
- stream_name, (long)htcop->pagesize);
- res = RES_BAD_ARG;
- goto error;
- }
-
- if(!IS_ALIGNED(htcop->pagesize, htcop->pagesize_os)) {
- log_err(htcop,
-"%s: invalid page size `%li'. The page size attribute must be aligned on the "
-"page size of the operating system , i.e. %lu.\n",
- stream_name, (long)htcop->pagesize, (unsigned long)htcop->pagesize_os);
- res = RES_BAD_ARG;
- goto error;
- }
-
- READ(&htcop->irregular_z, 1, "'irregular Z' flag");
- READ(htcop->definition, 4, "spatial and time definitions");
- if(htcop->definition[0] <= 0
- || htcop->definition[1] <= 0
- || htcop->definition[2] <= 0
- || htcop->definition[3] <= 0) {
- log_err(htcop,
-"%s: the spatial/time definition cannot be negative or null -- spatial "
-"definition: %i %i %i; time definition: %i\n",
- stream_name, SPLIT3(htcop->definition), htcop->definition[TIME]);
- res = RES_BAD_ARG;
- goto error;
- }
-
- READ(htcop->lower, 3, "lower position");
- READ(htcop->vxsz, 2, "XY voxel size ");
-
- nz = htcop->irregular_z ? (size_t)htcop->definition[Z] : 1;
- res = darray_double_resize(&htcop->vxsz_z, nz);
- if(res != RES_OK) {
- log_err(htcop,
- "%s: could not allocate memory to store the size of the voxels in Z.\n",
- stream_name);
- res = RES_MEM_ERR;
- goto error;
- }
- READ(darray_double_data_get(&htcop->vxsz_z), nz, "Z voxel size(s)");
- #undef READ
-
- map_len =
- (size_t)htcop->definition[X]
- * (size_t)htcop->definition[Y]
- * (size_t)htcop->definition[Z]
- * (size_t)htcop->definition[TIME]
- * sizeof(double);
- /* Align sizeof mapped data onto page size */
- map_len = ALIGN_SIZE(map_len, (size_t)htcop->pagesize);
-
- /* Ensure that offset is align on page size */
- offset = ALIGN_SIZE(ftell(stream), htcop->pagesize);
-
- #define MMAP(Var) { \
- htcop->Var = mmap(NULL, map_len, PROT_READ, MAP_PRIVATE|MAP_POPULATE, \
- fileno(stream), offset); \
- if(htcop->Var == MAP_FAILED) { \
- log_err(htcop, "%s: could not map the "STR(Var)" data -- %s.\n", \
- stream_name, strerror(errno)); \
- res = RES_IO_ERR; \
- goto error; \
- } \
- htcop->Var##_length = map_len; \
- } (void)0
- MMAP(RVT); offset += (off_t)map_len;
- MMAP(RCT); offset += (off_t)map_len;
- #undef MMAP
-
- CHK(fseek(stream, offset, SEEK_CUR) != -1);
-
-exit:
- return res;
-error:
- goto exit;
-}
-
-static void
-release_htcop(ref_T* ref)
-{
- struct htcop* htcop;
- ASSERT(ref);
- htcop = CONTAINER_OF(ref, struct htcop, ref);
- reset_htcop(htcop);
- darray_double_release(&htcop->vxsz_z);
- MEM_RM(htcop->allocator, htcop);
-}
-
-/*******************************************************************************
- * Exported functions
- ******************************************************************************/
-res_T
-htcop_create
- (struct logger* log,
- struct mem_allocator* mem_allocator,
- const int verbose,
- struct htcop** out_htcop)
-{
- struct htcop* htcop = NULL;
- struct mem_allocator* allocator = NULL;
- struct logger* logger = NULL;
- res_T res = RES_OK;
-
- if(!out_htcop) {
- res = RES_BAD_ARG;
- goto error;
- }
-
- allocator = mem_allocator ? mem_allocator : &mem_default_allocator;
- logger = log ? log : LOGGER_DEFAULT;
-
- htcop = MEM_CALLOC(allocator, 1, sizeof(struct htcop));
- if(!htcop) {
- if(verbose) {
- logger_print(logger, LOG_ERROR,
- "%s: could not allocate the HTCOP handler.\n", FUNC_NAME);
- res = RES_MEM_ERR;
- goto error;
- }
- }
- ref_init(&htcop->ref);
- htcop->allocator = allocator;
- htcop->logger = logger;
- htcop->verbose = verbose;
- htcop->pagesize_os = (size_t)sysconf(_SC_PAGESIZE);
- darray_double_init(htcop->allocator, &htcop->vxsz_z);
-
-exit:
- if(out_htcop) *out_htcop = htcop;
- return res;
-error:
- if(htcop) {
- HTCOP(ref_put(htcop));
- htcop = NULL;
- }
- goto exit;
-}
-
-res_T
-htcop_ref_get(struct htcop* htcop)
-{
- if(!htcop) return RES_BAD_ARG;
- ref_get(&htcop->ref);
- return RES_OK;
-}
-
-res_T
-htcop_ref_put(struct htcop* htcop)
-{
- if(!htcop) return RES_BAD_ARG;
- ref_put(&htcop->ref, release_htcop);
- return RES_OK;
-}
-
-res_T
-htcop_load(struct htcop* htcop, const char* path)
-{
- FILE* file = NULL;
- res_T res = RES_OK;
-
- if(!htcop || !path) {
- res = RES_BAD_ARG;
- goto error;
- }
-
- file = fopen(path, "r");
- if(!file) {
- log_err(htcop, "%s: error opening file `%S'.\n", FUNC_NAME, path);
- res = RES_IO_ERR;
- goto error;
- }
-
- res = load_stream(htcop, file, path);
- if(res != RES_OK) goto error;
-
-exit:
- if(file) fclose(file);
- return res;
-error:
- goto exit;
-}
-
-res_T
-htcop_load_stream(struct htcop* htcop, FILE* stream)
-{
- if(!htcop || !stream) return RES_BAD_ARG;
- return load_stream(htcop, stream, "<stream>");
-}
-
-res_T
-htcop_get_desc(const struct htcop* htcop, struct htcop_desc* desc)
-{
- if(!htcop || !desc) return RES_BAD_ARG;
- if(!htcop->RVT) return RES_BAD_ARG;
- desc->pagesize = (size_t)htcop->pagesize;
- desc->irregular_z = htcop->irregular_z != 0;
- desc->spatial_definition[0] = (size_t)htcop->definition[X];
- desc->spatial_definition[1] = (size_t)htcop->definition[Y];
- desc->spatial_definition[2] = (size_t)htcop->definition[Z];
- desc->time_definition = (size_t)htcop->definition[TIME];
- desc->lower[0] = htcop->lower[0];
- desc->lower[1] = htcop->lower[1];
- desc->lower[2] = htcop->lower[2];
- desc->vxsz_x = htcop->vxsz[0];
- desc->vxsz_y = htcop->vxsz[1];
- desc->vxsz_z = darray_double_cdata_get(&htcop->vxsz_z);
- desc->RCT = htcop->RCT;
- desc->RVT = htcop->RVT;
- return RES_OK;
-}
-
diff --git a/src/htcop.h b/src/htcop.h
@@ -1,139 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#ifndef HTCOP_H
-#define HTCOP_H
-
-#include <rsys/rsys.h>
-
-/* Library symbol management */
-#if defined(HTCOP_SHARED_BUILD) /* Build shared library */
- #define HTCOP_API extern EXPORT_SYM
-#elif defined(HTCOP_STATIC) /* Use/build static library */
- #define HTCOP_API extern LOCAL_SYM
-#else /* Use shared library */
- #define HTCOP_API extern IMPORT_SYM
-#endif
-
-/* Helper macro that asserts if the invocation of the htcop function `Func'
- * returns an error. One should use this macro on htcop function calls for
- * which no explicit error checking is performed */
-#ifndef NDEBUG
- #define HTCOP(Func) ASSERT(htcop_ ## Func == RES_OK)
-#else
- #define HTCOP(Func) htcop_ ## Func
-#endif
-
-/* Forward declaration of external data types */
-struct logger;
-struct mem_allocator;
-
-/* Forward declaration of opaque data types */
-struct htcop;
-
-struct htcop_desc {
- size_t pagesize; /* In bytes */
- int irregular_z;
-
- /* #coordinates along the X, Y and Z axis. X means for the West-Est axis, Y
- * is the South-North axis and finally Z is the altitude */
- size_t spatial_definition[3];
- size_t time_definition; /* Definition of the time */
-
- double lower[3]; /* Lower position of the grid */
- double vxsz_x; /* Voxel size in X */
- double vxsz_y; /* Voxel size in Y */
- const double* vxsz_z; /* Voxel size along Z */
-
- const double* RCT;
- const double* RVT;
-};
-#define HTCOP_DESC_NULL__ {0,-1,{0,0,0},0,{-1,-1,-1},-1,-1,NULL,NULL,NULL}
-static const struct htcop_desc HTCOP_DESC_NULL = HTCOP_DESC_NULL__;
-
-BEGIN_DECLS
-
-/*******************************************************************************
- * HTCOP API
- ******************************************************************************/
-HTCOP_API res_T
-htcop_create
-(struct logger* logger, /* NULL <=> use default logger */
- struct mem_allocator* allocator, /* NULL <=> use default allocator */
- const int verbose, /* Verbosity level */
- struct htcop** htcop);
-
-HTCOP_API res_T
-htcop_load
- (struct htcop* htcop,
- const char* path);
-
-HTCOP_API res_T
-htcop_load_stream
- (struct htcop* htcop,
- FILE* stream);
-
-HTCOP_API res_T
-htcop_ref_get
- (struct htcop* htcop);
-
-HTCOP_API res_T
-htcop_ref_put
- (struct htcop* htcop);
-
-HTCOP_API res_T
-htcop_get_desc
- (const struct htcop* htcop,
- struct htcop_desc* desc);
-
-/* Internal function */
-static FINLINE double
-htcop_dblgrid4D_at__
- (const double* grid,
- const struct htcop_desc* desc,
- size_t x, size_t y, size_t z, size_t t)
-{
- size_t row, slice, array;
- ASSERT(desc && grid);
- ASSERT(x < desc->spatial_definition[0]);
- ASSERT(y < desc->spatial_definition[1]);
- ASSERT(z < desc->spatial_definition[2]);
- ASSERT(t < desc->time_definition);
- row = desc->spatial_definition[0];
- slice = desc->spatial_definition[1] * row;
- array = desc->spatial_definition[2] * slice;
- return grid[t*array + z*slice + y*row + x];
-}
-
-static FINLINE double
-htcop_desc_RCT_at
- (const struct htcop_desc* desc,
- size_t x, size_t y, size_t z, size_t t)
-{
- return htcop_dblgrid4D_at__(desc->RCT, desc, x, y, z, t);
-}
-
-
-static FINLINE double
-htcop_desc_RVT_at
- (const struct htcop_desc* desc,
- size_t x, size_t y, size_t z, size_t t)
-{
- return htcop_dblgrid4D_at__(desc->RVT, desc, x, y, z, t);
-}
-
-END_DECLS
-
-#endif /* HTCOP_H */
diff --git a/src/htcp.c b/src/htcp.c
@@ -0,0 +1,335 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#define _POSIX_C_SOURCE 200809L /* mmap support */
+#define _DEFAULT_SOURCE 1 /* MAP_POPULATE support */
+
+#include "htcp.h"
+
+#include <rsys/dynamic_array_double.h>
+#include <rsys/logger.h>
+#include <rsys/ref_count.h>
+#include <rsys/mem_allocator.h>
+
+#include <errno.h>
+#include <string.h>
+#include <sys/mman.h>
+#include <unistd.h>
+
+enum { X, Y, Z, TIME }; /* Helper constants */
+
+struct htcp {
+ int64_t pagesize;
+ int8_t irregular_z;
+ int32_t definition[4];
+ double lower[3];
+ double vxsz[2]; /* Size of the voxels in X and Y */
+ struct darray_double vxsz_z; /* Size of the voxels along the Z dimension */
+
+ double* RCT; /* Mapped memory */
+ double* RVT; /* Mapped memory */
+ size_t RCT_length; /* In bytes */
+ size_t RVT_length; /* In bytes */
+
+ size_t pagesize_os; /* Page size of the os */
+ int verbose; /* Verbosity level */
+ struct logger* logger;
+ struct mem_allocator* allocator;
+ ref_T ref;
+};
+
+/*******************************************************************************
+ * Local functions
+ ******************************************************************************/
+static void
+log_msg
+ (const struct htcp* htcp,
+ const enum log_type stream,
+ const char* msg,
+ va_list vargs)
+{
+ ASSERT(htcp && msg);
+ if(htcp->verbose) {
+ res_T res; (void)res;
+ res = logger_vprint(htcp->logger, stream, msg, vargs);
+ ASSERT(res == RES_OK);
+ }
+}
+
+static void
+log_err(const struct htcp* htcp, const char* msg, ...)
+{
+ va_list vargs_list;
+ ASSERT(htcp && msg);
+ va_start(vargs_list, msg);
+ log_msg(htcp, LOG_ERROR, msg, vargs_list);
+ va_end(vargs_list);
+}
+
+static void
+reset_htcp(struct htcp* htcp)
+{
+ htcp->pagesize = 0;
+ htcp->irregular_z = 0;
+ htcp->definition[0] = 0;
+ htcp->definition[1] = 0;
+ htcp->definition[2] = 0;
+ htcp->lower[0] = -1;
+ htcp->lower[1] = -1;
+ htcp->lower[2] = -1;
+ htcp->vxsz[0] = -1;
+ htcp->vxsz[1] = -1;
+ darray_double_clear(&htcp->vxsz_z);
+ if(htcp->RCT) munmap(htcp->RCT, htcp->RCT_length);
+ if(htcp->RVT) munmap(htcp->RVT, htcp->RVT_length);
+ htcp->RCT = NULL;
+ htcp->RVT = NULL;
+ htcp->RCT_length = 0;
+ htcp->RVT_length = 0;
+}
+
+static res_T
+load_stream(struct htcp* htcp, FILE* stream, const char* stream_name)
+{
+ size_t nz = 0;
+ size_t map_len = 0;
+ off_t offset = 0;
+ res_T res = RES_OK;
+ ASSERT(htcp && stream && stream_name);
+
+ reset_htcp(htcp);
+
+ #define READ(Var, N, Name) { \
+ if(fread((Var), sizeof(*(Var)), (N), stream) != (N)) { \
+ fprintf(stderr, "%s: could not read the %s\n", stream_name, Name); \
+ res = RES_IO_ERR; \
+ goto error; \
+ } \
+ } (void)0
+ READ(&htcp->pagesize, 1, "page size");
+ if(!IS_POW2(htcp->pagesize)) {
+ log_err(htcp, "%s: invalid page size `%li'. It must be a power of two.\n",
+ stream_name, (long)htcp->pagesize);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ if(!IS_ALIGNED(htcp->pagesize, htcp->pagesize_os)) {
+ log_err(htcp,
+"%s: invalid page size `%li'. The page size attribute must be aligned on the "
+"page size of the operating system , i.e. %lu.\n",
+ stream_name, (long)htcp->pagesize, (unsigned long)htcp->pagesize_os);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ READ(&htcp->irregular_z, 1, "'irregular Z' flag");
+ READ(htcp->definition, 4, "spatial and time definitions");
+ if(htcp->definition[0] <= 0
+ || htcp->definition[1] <= 0
+ || htcp->definition[2] <= 0
+ || htcp->definition[3] <= 0) {
+ log_err(htcp,
+"%s: the spatial/time definition cannot be negative or null -- spatial "
+"definition: %i %i %i; time definition: %i\n",
+ stream_name, SPLIT3(htcp->definition), htcp->definition[TIME]);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ READ(htcp->lower, 3, "lower position");
+ READ(htcp->vxsz, 2, "XY voxel size ");
+
+ nz = htcp->irregular_z ? (size_t)htcp->definition[Z] : 1;
+ res = darray_double_resize(&htcp->vxsz_z, nz);
+ if(res != RES_OK) {
+ log_err(htcp,
+ "%s: could not allocate memory to store the size of the voxels in Z.\n",
+ stream_name);
+ res = RES_MEM_ERR;
+ goto error;
+ }
+ READ(darray_double_data_get(&htcp->vxsz_z), nz, "Z voxel size(s)");
+ #undef READ
+
+ map_len =
+ (size_t)htcp->definition[X]
+ * (size_t)htcp->definition[Y]
+ * (size_t)htcp->definition[Z]
+ * (size_t)htcp->definition[TIME]
+ * sizeof(double);
+ /* Align sizeof mapped data onto page size */
+ map_len = ALIGN_SIZE(map_len, (size_t)htcp->pagesize);
+
+ /* Ensure that offset is align on page size */
+ offset = ALIGN_SIZE(ftell(stream), htcp->pagesize);
+
+ #define MMAP(Var) { \
+ htcp->Var = mmap(NULL, map_len, PROT_READ, MAP_PRIVATE|MAP_POPULATE, \
+ fileno(stream), offset); \
+ if(htcp->Var == MAP_FAILED) { \
+ log_err(htcp, "%s: could not map the "STR(Var)" data -- %s.\n", \
+ stream_name, strerror(errno)); \
+ res = RES_IO_ERR; \
+ goto error; \
+ } \
+ htcp->Var##_length = map_len; \
+ } (void)0
+ MMAP(RVT); offset += (off_t)map_len;
+ MMAP(RCT); offset += (off_t)map_len;
+ #undef MMAP
+
+ CHK(fseek(stream, offset, SEEK_CUR) != -1);
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
+static void
+release_htcp(ref_T* ref)
+{
+ struct htcp* htcp;
+ ASSERT(ref);
+ htcp = CONTAINER_OF(ref, struct htcp, ref);
+ reset_htcp(htcp);
+ darray_double_release(&htcp->vxsz_z);
+ MEM_RM(htcp->allocator, htcp);
+}
+
+/*******************************************************************************
+ * Exported functions
+ ******************************************************************************/
+res_T
+htcp_create
+ (struct logger* log,
+ struct mem_allocator* mem_allocator,
+ const int verbose,
+ struct htcp** out_htcp)
+{
+ struct htcp* htcp = NULL;
+ struct mem_allocator* allocator = NULL;
+ struct logger* logger = NULL;
+ res_T res = RES_OK;
+
+ if(!out_htcp) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ allocator = mem_allocator ? mem_allocator : &mem_default_allocator;
+ logger = log ? log : LOGGER_DEFAULT;
+
+ htcp = MEM_CALLOC(allocator, 1, sizeof(struct htcp));
+ if(!htcp) {
+ if(verbose) {
+ logger_print(logger, LOG_ERROR,
+ "%s: could not allocate the HTCP handler.\n", FUNC_NAME);
+ res = RES_MEM_ERR;
+ goto error;
+ }
+ }
+ ref_init(&htcp->ref);
+ htcp->allocator = allocator;
+ htcp->logger = logger;
+ htcp->verbose = verbose;
+ htcp->pagesize_os = (size_t)sysconf(_SC_PAGESIZE);
+ darray_double_init(htcp->allocator, &htcp->vxsz_z);
+
+exit:
+ if(out_htcp) *out_htcp = htcp;
+ return res;
+error:
+ if(htcp) {
+ HTCP(ref_put(htcp));
+ htcp = NULL;
+ }
+ goto exit;
+}
+
+res_T
+htcp_ref_get(struct htcp* htcp)
+{
+ if(!htcp) return RES_BAD_ARG;
+ ref_get(&htcp->ref);
+ return RES_OK;
+}
+
+res_T
+htcp_ref_put(struct htcp* htcp)
+{
+ if(!htcp) return RES_BAD_ARG;
+ ref_put(&htcp->ref, release_htcp);
+ return RES_OK;
+}
+
+res_T
+htcp_load(struct htcp* htcp, const char* path)
+{
+ FILE* file = NULL;
+ res_T res = RES_OK;
+
+ if(!htcp || !path) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ file = fopen(path, "r");
+ if(!file) {
+ log_err(htcp, "%s: error opening file `%S'.\n", FUNC_NAME, path);
+ res = RES_IO_ERR;
+ goto error;
+ }
+
+ res = load_stream(htcp, file, path);
+ if(res != RES_OK) goto error;
+
+exit:
+ if(file) fclose(file);
+ return res;
+error:
+ goto exit;
+}
+
+res_T
+htcp_load_stream(struct htcp* htcp, FILE* stream)
+{
+ if(!htcp || !stream) return RES_BAD_ARG;
+ return load_stream(htcp, stream, "<stream>");
+}
+
+res_T
+htcp_get_desc(const struct htcp* htcp, struct htcp_desc* desc)
+{
+ if(!htcp || !desc) return RES_BAD_ARG;
+ if(!htcp->RVT) return RES_BAD_ARG;
+ desc->pagesize = (size_t)htcp->pagesize;
+ desc->irregular_z = htcp->irregular_z != 0;
+ desc->spatial_definition[0] = (size_t)htcp->definition[X];
+ desc->spatial_definition[1] = (size_t)htcp->definition[Y];
+ desc->spatial_definition[2] = (size_t)htcp->definition[Z];
+ desc->time_definition = (size_t)htcp->definition[TIME];
+ desc->lower[0] = htcp->lower[0];
+ desc->lower[1] = htcp->lower[1];
+ desc->lower[2] = htcp->lower[2];
+ desc->vxsz_x = htcp->vxsz[0];
+ desc->vxsz_y = htcp->vxsz[1];
+ desc->vxsz_z = darray_double_cdata_get(&htcp->vxsz_z);
+ desc->RCT = htcp->RCT;
+ desc->RVT = htcp->RVT;
+ return RES_OK;
+}
+
diff --git a/src/htcp.h b/src/htcp.h
@@ -0,0 +1,139 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#ifndef HTCP_H
+#define HTCP_H
+
+#include <rsys/rsys.h>
+
+/* Library symbol management */
+#if defined(HTCP_SHARED_BUILD) /* Build shared library */
+ #define HTCP_API extern EXPORT_SYM
+#elif defined(HTCP_STATIC) /* Use/build static library */
+ #define HTCP_API extern LOCAL_SYM
+#else /* Use shared library */
+ #define HTCP_API extern IMPORT_SYM
+#endif
+
+/* Helper macro that asserts if the invocation of the htcp function `Func'
+ * returns an error. One should use this macro on htcp function calls for
+ * which no explicit error checking is performed */
+#ifndef NDEBUG
+ #define HTCP(Func) ASSERT(htcp_ ## Func == RES_OK)
+#else
+ #define HTCP(Func) htcp_ ## Func
+#endif
+
+/* Forward declaration of external data types */
+struct logger;
+struct mem_allocator;
+
+/* Forward declaration of opaque data types */
+struct htcp;
+
+struct htcp_desc {
+ size_t pagesize; /* In bytes */
+ int irregular_z;
+
+ /* #coordinates along the X, Y and Z axis. X means for the West-Est axis, Y
+ * is the South-North axis and finally Z is the altitude */
+ size_t spatial_definition[3];
+ size_t time_definition; /* Definition of the time */
+
+ double lower[3]; /* Lower position of the grid */
+ double vxsz_x; /* Voxel size in X */
+ double vxsz_y; /* Voxel size in Y */
+ const double* vxsz_z; /* Voxel size along Z */
+
+ const double* RCT;
+ const double* RVT;
+};
+#define HTCP_DESC_NULL__ {0,-1,{0,0,0},0,{-1,-1,-1},-1,-1,NULL,NULL,NULL}
+static const struct htcp_desc HTCP_DESC_NULL = HTCP_DESC_NULL__;
+
+BEGIN_DECLS
+
+/*******************************************************************************
+ * HTCP API
+ ******************************************************************************/
+HTCP_API res_T
+htcp_create
+(struct logger* logger, /* NULL <=> use default logger */
+ struct mem_allocator* allocator, /* NULL <=> use default allocator */
+ const int verbose, /* Verbosity level */
+ struct htcp** htcp);
+
+HTCP_API res_T
+htcp_load
+ (struct htcp* htcp,
+ const char* path);
+
+HTCP_API res_T
+htcp_load_stream
+ (struct htcp* htcp,
+ FILE* stream);
+
+HTCP_API res_T
+htcp_ref_get
+ (struct htcp* htcp);
+
+HTCP_API res_T
+htcp_ref_put
+ (struct htcp* htcp);
+
+HTCP_API res_T
+htcp_get_desc
+ (const struct htcp* htcp,
+ struct htcp_desc* desc);
+
+/* Internal function */
+static FINLINE double
+htcp_dblgrid4D_at__
+ (const double* grid,
+ const struct htcp_desc* desc,
+ size_t x, size_t y, size_t z, size_t t)
+{
+ size_t row, slice, array;
+ ASSERT(desc && grid);
+ ASSERT(x < desc->spatial_definition[0]);
+ ASSERT(y < desc->spatial_definition[1]);
+ ASSERT(z < desc->spatial_definition[2]);
+ ASSERT(t < desc->time_definition);
+ row = desc->spatial_definition[0];
+ slice = desc->spatial_definition[1] * row;
+ array = desc->spatial_definition[2] * slice;
+ return grid[t*array + z*slice + y*row + x];
+}
+
+static FINLINE double
+htcp_desc_RCT_at
+ (const struct htcp_desc* desc,
+ size_t x, size_t y, size_t z, size_t t)
+{
+ return htcp_dblgrid4D_at__(desc->RCT, desc, x, y, z, t);
+}
+
+
+static FINLINE double
+htcp_desc_RVT_at
+ (const struct htcp_desc* desc,
+ size_t x, size_t y, size_t z, size_t t)
+{
+ return htcp_dblgrid4D_at__(desc->RVT, desc, x, y, z, t);
+}
+
+END_DECLS
+
+#endif /* HTCP_H */
diff --git a/src/les2htcop.c b/src/les2htcop.c
@@ -1,790 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#define _POSIX_C_SOURCE 200112L /* getopt/close support */
-
-#include "les2htcop.h"
-
-#include <rsys/cstr.h>
-#include <rsys/math.h>
-#include <rsys/mem_allocator.h>
-#include <rsys/rsys.h>
-
-#include <alloca.h>
-#include <netcdf.h>
-#include <string.h>
-#include <unistd.h> /* getopt */
-#include <fcntl.h> /* open */
-#include <sys/stat.h> /* S_IRUSR & S_IWUSR */
-
-/*******************************************************************************
- * Command arguments
- ******************************************************************************/
-struct args {
- const char* output;
- const char* input;
- long pagesize;
- int force_overwrite;
- int check;
- int no_output;
- int quit; /* Quit the application */
-};
-#define ARGS_DEFAULT__ {NULL,NULL,4096,0,0,0,0}
-static const struct args ARGS_DEFAULT = ARGS_DEFAULT__;
-
-static void
-print_help(const char* cmd)
-{
- ASSERT(cmd);
-
- printf(
-"Usage: %s -i INPUT [OPTIONS]\n"
-"Convert the LES data stored into INPUT from NetCDF to the htcop fileformat.\n\n",
- cmd);
- printf(
-" -c advanced check of the validity of the submitted LES file\n"
-" with respect to the converter pre-requisites. Note that\n"
-" this option can increase significantly the conversion\n"
-" time.\n");
- printf(
-" -f overwrite the OUTPUT file if it already exists.\n");
- printf(
-" -h display this help and exit.\n");
- printf(
-" -i INPUT path of the LES file to convert.\n");
- printf(
-" -o OUTPUT write results to OUTPUT. If not defined, write results to\n"
-" standard output.\n");
- printf(
-" -p targeted page size in bytes. Must be a power of 2. The size\n"
-" of the converted LES data and their address into output are\n"
-" aligned according to this size. By default, the page size\n"
-" is %li Bytes.\n",
- ARGS_DEFAULT.pagesize);
- printf(
-" -q do not write results to OUTPUT.\n");
- printf(
-" -v display version information and exit.\n");
- printf("\n");
- printf(
-"%s (C) 2018 |Meso|Star>. This is free software released under the GNU GPL\n"
-"license, version 3 or later. You are free to change or redistribute it under\n"
-"certain conditions <http://gnu.org/licenses/gpl.html>.\n", cmd);
-}
-
-static void
-args_release(struct args* args)
-{
- ASSERT(args);
- *args = ARGS_DEFAULT;
-}
-
-static res_T
-args_init(struct args* args, const int argc, char** argv)
-{
- int opt;
- res_T res = RES_OK;
- ASSERT(args && argc && argv);
-
- while((opt = getopt(argc, argv, "cfhi:o:p:qv")) != -1) {
- switch(opt) {
- case 'c': args->check = 1; break;
- case 'f': args->force_overwrite = 1; break;
- case 'h':
- print_help(argv[0]);
- args_release(args);
- args->quit = 1;
- goto exit;
- case 'i': args->input = optarg; break;
- case 'o': args->output = optarg; break;
- case 'p':
- res = cstr_to_long(optarg, &args->pagesize);
- if(res == RES_OK && !IS_POW2(args->pagesize)) res = RES_BAD_ARG;
- break;
- case 'q': args->no_output = 1; break;
- case 'v':
- printf("%s %d.%d.%d\n",
- argv[0],
- LES2HTLES_VERSION_MAJOR,
- LES2HTLES_VERSION_MINOR,
- LES2HTLES_VERSION_PATCH);
- args->quit = 1;
- goto exit;
- default: RES_BAD_ARG; break;
- }
- if(res != RES_OK) {
- if(optarg) {
- fprintf(stderr, "%s: invalid option argumet '%s' -- '%c'\n",
- argv[0], optarg, opt);
- }
- goto error;
- }
- }
-
- if(!args->input) {
- fprintf(stderr, "Missing input file.\n");
- res = RES_BAD_ARG;
- goto error;
- }
-
- if(args->no_output) args->output = NULL;
-
-exit:
- return res;
-error:
- args_release(args);
- goto exit;
-}
-
-/*******************************************************************************
- * Grid header
- ******************************************************************************/
-struct grid {
- int32_t nx, ny, nz, ntimes; /* Definition */
- int8_t is_z_irregular;
- double lower[3]; /* Lower position of the grid */
- double vxsz_x; /* Size of the voxel is X */
- double vxsz_y; /* Size of the voxel in Y */
- double* vxsz_z; /* Size of the voxel in Z */
-};
-
-#define GRID_NULL__ {0,0,0,0,0,{0,0,0},0,0,NULL}
-static const struct grid GRID_NULL = GRID_NULL__;
-
-static void
-grid_release(struct grid* grid)
-{
- CHK(grid);
- if(grid->vxsz_z) mem_rm(grid->vxsz_z);
-}
-
-/*******************************************************************************
- * NetCDF helper functions and macros
- ******************************************************************************/
-#define INVALID_ID -1
-
-#define NC(Func) { \
- const int err__ = nc_ ## Func; \
- if(err__ != NC_NOERR) { \
- fprintf(stderr, "error:%i:%s\n", __LINE__, ncerr_to_str(err__)); \
- abort(); \
- } \
- } (void)0
-
-static INLINE const char*
-ncerr_to_str(const int ncerr)
-{
- const char* str = "NC_ERR_<UNKNOWN>";
- switch(ncerr) {
- case NC_EBADGRPID: str = "NC_EBADGRPID"; break;
- case NC_EBADID: str = "NC_EBADID"; break;
- case NC_EBADNAME: str = "NC_EBADNAME"; break;
- case NC_ECHAR: str = "NC_ECHAR"; break;
- case NC_EDIMMETA: str = "NC_EDIMMETA"; break;
- case NC_EHDFERR: str = "NC_EHDFERR"; break;
- case NC_ENOMEM: str = "NC_ENOMEM"; break;
- case NC_ENOTATT: str = "NC_ENOTATT"; break;
- case NC_ENOTVAR: str = "NC_ENOTVAR"; break;
- case NC_ERANGE: str = "NC_ERANGE"; break;
- case NC_NOERR: str = "NC_NOERR"; break;
- }
- return str;
-}
-
-static INLINE const char*
-nctype_to_str(const nc_type type)
-{
- const char* str = "NC_TYPE_<UNKNOWN>";
- switch(type) {
- case NC_NAT: str = "NC_NAT"; break;
- case NC_BYTE: str = "NC_BYTE"; break;
- case NC_CHAR: str = "NC_CHAR"; break;
- case NC_SHORT: str = "NC_SHORT"; break;
- case NC_LONG: str = "NC_LONG"; break;
- case NC_FLOAT: str = "NC_FLOAT"; break;
- case NC_DOUBLE: str = "NC_DOUBLE"; break;
- case NC_UBYTE: str = "NC_UBYTE"; break;
- case NC_USHORT: str = "NC_USHORT"; break;
- case NC_UINT: str = "NC_UINT"; break;
- case NC_INT64: str = "NC_INT64"; break;
- case NC_UINT64: str = "NC_UINT64"; break;
- case NC_STRING: str = "NC_STRING"; break;
- default: FATAL("Unreachable code.\n"); break;
- }
- return str;
-}
-
-static INLINE size_t
-sizeof_nctype(const nc_type type)
-{
- size_t sz;
- switch(type) {
- case NC_BYTE:
- case NC_CHAR:
- case NC_UBYTE:
- sz = 1; break;
- case NC_SHORT:
- case NC_USHORT:
- sz = 2; break;
- case NC_FLOAT:
- case NC_INT:
- case NC_UINT:
- sz = 4; break;
- case NC_DOUBLE:
- case NC_INT64:
- case NC_UINT64:
- sz = 8; break;
- default: FATAL("Unreachable cde.\n"); break;
- }
- return sz;
-}
-
-/*******************************************************************************
- * Helper functions
- ******************************************************************************/
-static res_T
-open_output_stream(const char* path, const int force_overwrite, FILE** stream)
-{
- int fd = -1;
- FILE* fp = NULL;
- res_T res = RES_OK;
- ASSERT(path);
-
- if(force_overwrite) {
- fp = fopen(path, "w");
- if(!fp) {
- fprintf(stderr, "Could not open the output file '%s'.\n", path);
- goto error;
- }
- } else {
- fd = open(path, O_CREAT|O_WRONLY|O_EXCL|O_TRUNC, S_IRUSR|S_IWUSR);
- if(fd >= 0) {
- fp = fdopen(fd, "w");
- if(fp == NULL) {
- fprintf(stderr, "Could not open the output file '%s'.\n", path);
- goto error;
- }
- } else if(errno == EEXIST) {
- fprintf(stderr,
- "The output file '%s' already exists. Use -f to overwrite it.\n", path);
- goto error;
- } else {
- fprintf(stderr,
- "Unexpected error while opening the output file `'%s'.\n", path);
- goto error;
- }
- }
-
-exit:
- *stream = fp;
- return res;
-error:
- res = RES_IO_ERR;
- if(fp) {
- CHK(fclose(fp) == 0);
- fp = NULL;
- } else if(fd >= 0) {
- CHK(close(fd) == 0);
- }
- goto exit;
-}
-
-
-static res_T
-setup_definition
- (int nc, int32_t* nx, int32_t* ny, int32_t* nz, int32_t* ntimes)
-{
- size_t len[4];
- int dimids[4];
- int id;
- int ndims;
- int err = NC_NOERR;
- res_T res = RES_OK;
- ASSERT(nx && ny && nz && ntimes);
-
- err = nc_inq_varid(nc, "RCT", &id);
- if(err != NC_NOERR) {
- fprintf(stderr, "Could not inquire the 'RCT' variable -- %s\n",
- ncerr_to_str(err));
- res = RES_BAD_ARG;
- goto error;
- }
-
- NC(inq_varndims(nc, id, &ndims));
- if(ndims != 4) {
- fprintf(stderr, "The dimension of the 'RCT' variable must be 4.\n");
- res = RES_BAD_ARG;
- goto error;
- }
-
- NC(inq_vardimid(nc, id, dimids));
- NC(inq_dimlen(nc, dimids[0], len+0));
- NC(inq_dimlen(nc, dimids[1], len+1));
- NC(inq_dimlen(nc, dimids[2], len+2));
- NC(inq_dimlen(nc, dimids[3], len+3));
-
- *nx = (int32_t)len[3];
- *ny = (int32_t)len[2];
- *nz = (int32_t)len[1];
- *ntimes = (int32_t)len[0];
-
- if(*nx < 1 || *ny < 1 || *nz < 1 || *ntimes < 1) {
- fprintf(stderr,
- "The spatial and time definitions cannot be null.\n"
- " #x = %i; #y = %i; #z = %i; #time = %i\n",
- *nx, *ny, *nz, *ntimes);
- res = RES_BAD_ARG;
- goto error;
- }
-
-exit:
- return res;
-error:
- goto exit;
-}
-
-static res_T
-setup_regular_dimension
- (const int nc,
- const char* var_name,
- double* voxel_size,
- double* lower_pos,
- const int check)
-{
- char* mem = NULL;
- size_t len, start, i;
- double vxsz;
- nc_type type;
- int varid, dimid, ndims;
- int err = NC_NOERR;
- res_T res = RES_OK;
- ASSERT(nc && var_name && voxel_size && lower_pos);
-
- err = nc_inq_varid(nc, var_name, &varid);
- if(err != NC_NOERR) {
- fprintf(stderr, "Could not inquire the '%s' variable -- %s\n",
- var_name, ncerr_to_str(err));
- res = RES_BAD_ARG;
- goto error;
- }
- NC(inq_varndims(nc, varid, &ndims));
- if(ndims != 1) {
- fprintf(stderr, "The dimension of the '%s' variable must be 1 .\n", var_name);
- res = RES_BAD_ARG;
- goto error;
- }
-
- NC(inq_vardimid(nc, varid, &dimid));
- NC(inq_dimlen(nc, dimid, &len));
- ASSERT(len > 0);
-
- NC(inq_vartype(nc, varid, &type));
- if(type != NC_DOUBLE && type != NC_FLOAT) {
- fprintf(stderr, "The type of the '%s' variable must be float or double.\n",
- var_name);
- res = RES_BAD_ARG;
- goto error;
- }
-
- #define AT(Mem, Id) (type == NC_FLOAT ? ((float*)Mem)[Id] : ((double*)Mem)[Id])
- if(!check) {
- char* raw = alloca(8);
- start = 0;
- len = 1;
- NC(get_vara(nc, varid, &start, &len, raw));
- vxsz = AT(raw, 0) * 2; /* Assume that the grid starts at 0 */
- } else {
- mem = mem_alloc(len*sizeof_nctype(type));
- if(!mem) {
- fprintf(stderr, "Could not allocate memory for the '%s' variable.\n", var_name);
- res = RES_MEM_ERR;
- goto error;
- }
- start = 0;
- NC(get_vara(nc, varid, &start, &len, mem));
- vxsz = AT(mem, 0) * 2; /* Assume that the grid starts from 0 */
- }
-
- if(vxsz <= 0) {
- fprintf(stderr, "The '%s' variable can't have negative or null values.\n",
- var_name);
- res = RES_BAD_ARG;
- goto error;
- }
-
- if(check) {
- /* Check that the dimension is regular */
- FOR_EACH(i, 1, len) {
- if(!eq_eps(AT(mem, i) - AT(mem, i-1), vxsz, 1.e-6)) {
- fprintf(stderr, "The voxel size of the '%s' variable must be regular.\n",
- var_name);
- res = RES_BAD_ARG;
- goto error;
- }
- }
- }
- #undef AT
-
- *voxel_size = vxsz;
- *lower_pos = 0;
-exit:
- if(mem) mem_rm(mem);
- return res;
-error:
- goto exit;
-}
-
-/* Return 1 if Z is irregular */
-static int
-setup_Z_dimension
- (const int nc,
- double** voxel_size,
- double* lower_pos,
- int8_t* is_z_irregular,
- const int check)
-{
- enum { Z, X, Y, NDIMS }; /* Helper constants */
- char* mem = NULL; /* Memory where NetCDF data are copied */
- char* mem2 = NULL; /* Temporary memory used to check the VLEV data */
- size_t len[NDIMS];
- size_t start[NDIMS];
- size_t dim_len[NDIMS];
- size_t z; /* Id over the Z dimension */
- size_t i; /* Common id */
- double vxsz; /* Voxel size */
- double* pvxsz = NULL; /* Pointer toward voxel sizes */
- nc_type type; /* Type of the NetCDF VLEV variable */
- int varid; /* NetCDF id of the VLEV variable */
- int dimids[NDIMS]; /* NetCDF id of the VLEV dimensions */
- int irregular; /* Boolean defining if the Z dimension is irregular */
- int ndims; /* #dims of the VLEV variable */
- int err = NC_NOERR;
- res_T res = RES_OK;
- ASSERT(nc && voxel_size && lower_pos && is_z_irregular);
-
- /* Retrieve the VLEV variable */
- err = nc_inq_varid(nc, "VLEV", &varid);
- if(err != NC_NOERR) {
- fprintf(stderr, "Could not inquire the 'VLEV' variable -- %s\n",
- ncerr_to_str(err));
- res = RES_BAD_ARG;
- goto error;
- }
-
- /* Inquire the VLEV dimensionality */
- NC(inq_varndims(nc, varid, &ndims));
- if(ndims != NDIMS) {
- fprintf(stderr, "The dimension of the 'VLEV' variable must be %i .\n",
- NDIMS);
- res = RES_BAD_ARG;
- goto error;
- }
- NC(inq_vardimid(nc, varid, dimids));
- FOR_EACH(i, 0, NDIMS) NC(inq_dimlen(nc, dimids[i], len+i));
-
- /* Inquire the type of the VLEV variable */
- NC(inq_vartype(nc, varid, &type));
- if(type != NC_DOUBLE && type != NC_FLOAT) {
- fprintf(stderr, "The type of 'VLEV' variable must be float or double.\n");
- res = RES_BAD_ARG;
- goto error;
- }
-
- /* Allocate te memory space where a Z column of the VLEV variable is going to
- * be copied */
- mem = mem_alloc(len[Z]*sizeof_nctype(type));
- if(!mem) {
- fprintf(stderr, "Could not allocate memory for the 'VLEV' variable.\n");
- res = RES_MEM_ERR;
- goto error;
- }
-
- /* Get only one Z column */
- start[X] = 0; dim_len[X] = 1;
- start[Y] = 0; dim_len[Y] = 1;
- start[Z] = 0; dim_len[Z] = len[Z];
- NC(get_vara(nc, varid, start, dim_len, mem));
-
- #define AT(Mem, Id) (type == NC_FLOAT ? ((float*)Mem)[Id] : ((double*)Mem)[Id])
-
- /* Assume that the Z dimension starts from 0 */
- vxsz = AT(mem, 0) * 2.0;
- if(vxsz <= 0) {
- fprintf(stderr, "The 'VLEV' variable can't have negative or null values.\n");
- res = RES_BAD_ARG;
- goto error;
- }
-
- /* Check if the Z dimension is regular */
- FOR_EACH(z, 1, len[Z]) {
- if(!eq_eps(AT(mem, z) - AT(mem, z-1), vxsz, 1.e-6)) break;
- }
- irregular = (z != len[Z]);
-
- /* Setup the size of the voxel */
- pvxsz = mem_alloc(sizeof(double) * (irregular ? len[Z] : 1));
- if(!pvxsz) {
- fprintf(stderr,
- "Could not allocate memory for the voxel size along the Z dimension.\n");
- res = RES_MEM_ERR;
- goto error;
- }
- pvxsz[0] = vxsz;
-
- if(irregular) {
- FOR_EACH(z, 1, len[Z]) {
- pvxsz[z] = (AT(mem, z) - (AT(mem, z-1) + pvxsz[z-1]*0.5)) * 2.0;
- if(pvxsz[z] <= 0) {
- fprintf(stderr,
- "The 'VLEV' variable can't have negative or null values.\n");
- res = RES_BAD_ARG;
- goto error;
- }
- }
- }
-
- /* Check that the others columns have the same voxel size */
- if(check) {
- size_t ix, iy;
-
- mem2 = mem_alloc(len[Z]*sizeof_nctype(type));
- if(!mem2) {
- fprintf(stderr,
- "Could not allocate memory to check the 'VLEV' variable.\n");
- res = RES_MEM_ERR;
- goto error;
- }
-
- FOR_EACH(iy, 0, len[Y]) {
- start[Y] = iy;
- FOR_EACH(ix, 1, len[X]) {
- start[X] = ix;
- NC(get_vara(nc, varid, start, dim_len, mem2));
- FOR_EACH(z, 0, len[Z]) {
- if(!eq_eps(AT(mem2, z), AT(mem, z), 1.e-6)) {
- fprintf(stderr,
- "The Z columns of the 'VLEV' variable must be the same.\n");
- res = RES_BAD_ARG;
- goto error;
- }
- }
- }
- }
- }
- #undef AT
-
- *lower_pos = 0;
- *voxel_size = pvxsz;
- *is_z_irregular = (int8_t)irregular;
-
-exit:
- if(mem2) mem_rm(mem2);
- if(mem) mem_rm(mem);
- return res;
-error:
- goto exit;
-}
-
-static res_T
-write_data(int nc, const char* var_name, FILE* stream)
-{
- enum { TIME, Z, X, Y, NDIMS }; /* Helper constants */
-
- char* mem = NULL; /* Memory where NetCDF data are copied */
- char* mem2 = NULL; /* Tmp memory where NetCDF data are converted */
- size_t start[NDIMS];
- size_t dim_len[NDIMS];
- size_t len[NDIMS];
- size_t slice_len; /* #elements per slice the 3D grid */
- size_t i; /* Common id */
- size_t t; /* Id over the time dimension */
- size_t z; /* Id over the Z dimension */
- nc_type type;
- int varid; /* NetCDF id of the data to write */
- int ndims; /* #dimensions of the NetCDF data to write */
- int dimids[NDIMS]; /* NetCDF id of the data dimension */
- int err = NC_NOERR;
- res_T res = RES_OK;
- CHK(var_name);
-
- /* Retrieve the NetCDF id of the variable */
- err = nc_inq_varid(nc, var_name, &varid);
- if(err != NC_NOERR) {
- fprintf(stderr, "Could not inquire the '%s' variable -- %s\n",
- var_name, ncerr_to_str(err));
- res = RES_BAD_ARG;
- goto error;
- }
-
- /* Inquire the dimensions of the variable */
- NC(inq_varndims(nc, varid, &ndims));
- if(ndims != NDIMS) {
- fprintf(stderr, "The dimension of the '%s' variable must be %i .\n",
- var_name, NDIMS);
- res = RES_BAD_ARG;
- goto error;
- }
- NC(inq_vardimid(nc, varid, dimids));
- FOR_EACH(i, 0, NDIMS) NC(inq_dimlen(nc, dimids[i], len+i));
-
- /* Inquire the type of the variable */
- NC(inq_vartype(nc, varid, &type));
- if(type != NC_DOUBLE && type != NC_FLOAT) {
- fprintf(stderr,
- "The type of the '%s' variable cannot be %s. Expecting floating point data.\n",
- var_name, nctype_to_str(type));
- res = RES_BAD_ARG;
- goto error;
- }
-
- /* Alloc the memory space where the variable data wille be copied. Note that,
- * in order to reduce memory consumption, only one slice is read at a given
- * time and Z position. */
- slice_len = len[X]*len[Y]; /* # elements per slice */
- mem = mem_alloc(slice_len*sizeof_nctype(type));
-
- /* Allocate a temporary memory space to store the variable data converted in
- * double precision if they are stored in single precision. */
- mem2 = type==NC_DOUBLE ? mem : mem_alloc(slice_len*sizeof_nctype(NC_DOUBLE));
-
- if(!mem || !mem2) {
- fprintf(stderr, "Could not allocate memory for the '%s' variable.\n",
- var_name);
- res = RES_MEM_ERR;
- goto error;
- }
-
- /* Read the whole slice */
- start[X] = 0; dim_len[X] = len[X];
- start[Y] = 0; dim_len[Y] = len[Y];
-
- FOR_EACH(t, 0, len[TIME]) {
- start[TIME] = t, dim_len[TIME] = 1;
- FOR_EACH(z, 0, len[Z]) {
- size_t n;
- start[Z] = z; dim_len[Z] = 1;
-
- NC(get_vara(nc, varid, start, dim_len, mem)); /* Fetch the slice data */
-
- /* In order to be compliant with the htgop fileformat, single precision
- * data must be converted in double precision */
- if(type == NC_FLOAT) {
- FOR_EACH(i, 0, slice_len) ((double*)mem2)[i] = ((float*)mem)[i];
- }
-
- /* Write data */
- n = fwrite(mem2, sizeof_nctype(NC_DOUBLE), slice_len, stream);
- if(n != slice_len) {
- fprintf(stderr, "Error writing data of the '%s' variable -- '%s'.\n",
- var_name, strerror(errno));
- res = RES_IO_ERR;
- goto error;
- }
- }
- }
-
-exit:
- if(mem2 != mem) mem_rm(mem2);
- if(mem) mem_rm(mem);
- return res;
-error:
- goto exit;
-}
-
-/*******************************************************************************
- * Program
- ******************************************************************************/
-int
-main(int argc, char** argv)
-{
- FILE* stream = stdout;
- struct grid grid = GRID_NULL;
- struct args args = ARGS_DEFAULT;
- int err = 0;
- int err_nc = NC_NOERR;
- int nc = INVALID_ID;
- res_T res = RES_OK;
-
- res = args_init(&args, argc, argv);
- if(res != RES_OK) goto error;
- if(args.quit) goto exit;
-
- if(args.output) {
- res = open_output_stream(args.output, args.force_overwrite, &stream);
- if(res != RES_OK) goto error;
- }
-
- err_nc = nc_open(args.input, NC_WRITE|NC_MMAP, &nc);
- if(err_nc != NC_NOERR) {
- fprintf(stderr, "Error opening file `%s' -- %s.\n",
- args.input, ncerr_to_str(err_nc));
- goto exit;
- }
-
- #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0
- CALL(setup_definition(nc, &grid.nx, &grid.ny, &grid.nz, &grid.ntimes));
- CALL(setup_regular_dimension
- (nc, "W_E_direction", &grid.vxsz_x, grid.lower+0, args.check));
- CALL(setup_regular_dimension
- (nc, "S_N_direction", &grid.vxsz_y, grid.lower+1, args.check));
- CALL(setup_Z_dimension
- (nc, &grid.vxsz_z, grid.lower+2, &grid.is_z_irregular, args.check));
- #undef CALL
-
- #define WRITE(Var, N, Name) { \
- if(fwrite((Var), sizeof(*(Var)), (N), stream) != (N)) { \
- fprintf(stderr, "Error writing the %s.\n", (Name)); \
- res = RES_IO_ERR; \
- goto error; \
- } \
- } (void)0
- if(!args.no_output) {
- const int64_t pagesize = (int64_t)args.pagesize;
-
- WRITE(&pagesize, 1, "page size");
- WRITE(&grid.is_z_irregular, 1, "'irregular' Z boolean");
- WRITE(&grid.nx, 1, "X definition");
- WRITE(&grid.ny, 1, "Y definition");
- WRITE(&grid.nz, 1, "Z definition");
- WRITE(&grid.ntimes, 1, "time definition");
- WRITE(grid.lower, 3, "lower position");
- WRITE(&grid.vxsz_x, 1, "X voxel size");
- WRITE(&grid.vxsz_y, 1, "Y voxel size");
- WRITE(grid.vxsz_z, grid.is_z_irregular?(size_t)grid.nz:1, "Z voxel size(s)");
- /* padding */
- fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
- write_data(nc, "RVT", stream);
- /* padding */
- fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
- write_data(nc, "RCT", stream);
- /*padding */
- fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
- }
- #undef WRITE
-
-exit:
- grid_release(&grid);
- if(nc != INVALID_ID) NC(close(nc));
- if(mem_allocated_size() != 0) {
- fprintf(stderr, "Memory leaks: %lu Bytes\n",
- (unsigned long)mem_allocated_size());
- err = -1;
- }
- return err;
-error:
- err = -1;
- goto exit;
-}
-
diff --git a/src/les2htcp.c b/src/les2htcp.c
@@ -0,0 +1,790 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#define _POSIX_C_SOURCE 200112L /* getopt/close support */
+
+#include "les2htcp.h"
+
+#include <rsys/cstr.h>
+#include <rsys/math.h>
+#include <rsys/mem_allocator.h>
+#include <rsys/rsys.h>
+
+#include <alloca.h>
+#include <netcdf.h>
+#include <string.h>
+#include <unistd.h> /* getopt */
+#include <fcntl.h> /* open */
+#include <sys/stat.h> /* S_IRUSR & S_IWUSR */
+
+/*******************************************************************************
+ * Command arguments
+ ******************************************************************************/
+struct args {
+ const char* output;
+ const char* input;
+ long pagesize;
+ int force_overwrite;
+ int check;
+ int no_output;
+ int quit; /* Quit the application */
+};
+#define ARGS_DEFAULT__ {NULL,NULL,4096,0,0,0,0}
+static const struct args ARGS_DEFAULT = ARGS_DEFAULT__;
+
+static void
+print_help(const char* cmd)
+{
+ ASSERT(cmd);
+
+ printf(
+"Usage: %s -i INPUT [OPTIONS]\n"
+"Convert the LES data stored into INPUT from NetCDF to the htcp fileformat.\n\n",
+ cmd);
+ printf(
+" -c advanced check of the validity of the submitted LES file\n"
+" with respect to the converter pre-requisites. Note that\n"
+" this option can increase significantly the conversion\n"
+" time.\n");
+ printf(
+" -f overwrite the OUTPUT file if it already exists.\n");
+ printf(
+" -h display this help and exit.\n");
+ printf(
+" -i INPUT path of the LES file to convert.\n");
+ printf(
+" -o OUTPUT write results to OUTPUT. If not defined, write results to\n"
+" standard output.\n");
+ printf(
+" -p targeted page size in bytes. Must be a power of 2. The size\n"
+" of the converted LES data and their address into output are\n"
+" aligned according to this size. By default, the page size\n"
+" is %li Bytes.\n",
+ ARGS_DEFAULT.pagesize);
+ printf(
+" -q do not write results to OUTPUT.\n");
+ printf(
+" -v display version information and exit.\n");
+ printf("\n");
+ printf(
+"%s (C) 2018 |Meso|Star>. This is free software released under the GNU GPL\n"
+"license, version 3 or later. You are free to change or redistribute it under\n"
+"certain conditions <http://gnu.org/licenses/gpl.html>.\n", cmd);
+}
+
+static void
+args_release(struct args* args)
+{
+ ASSERT(args);
+ *args = ARGS_DEFAULT;
+}
+
+static res_T
+args_init(struct args* args, const int argc, char** argv)
+{
+ int opt;
+ res_T res = RES_OK;
+ ASSERT(args && argc && argv);
+
+ while((opt = getopt(argc, argv, "cfhi:o:p:qv")) != -1) {
+ switch(opt) {
+ case 'c': args->check = 1; break;
+ case 'f': args->force_overwrite = 1; break;
+ case 'h':
+ print_help(argv[0]);
+ args_release(args);
+ args->quit = 1;
+ goto exit;
+ case 'i': args->input = optarg; break;
+ case 'o': args->output = optarg; break;
+ case 'p':
+ res = cstr_to_long(optarg, &args->pagesize);
+ if(res == RES_OK && !IS_POW2(args->pagesize)) res = RES_BAD_ARG;
+ break;
+ case 'q': args->no_output = 1; break;
+ case 'v':
+ printf("%s %d.%d.%d\n",
+ argv[0],
+ LES2HTLES_VERSION_MAJOR,
+ LES2HTLES_VERSION_MINOR,
+ LES2HTLES_VERSION_PATCH);
+ args->quit = 1;
+ goto exit;
+ default: RES_BAD_ARG; break;
+ }
+ if(res != RES_OK) {
+ if(optarg) {
+ fprintf(stderr, "%s: invalid option argumet '%s' -- '%c'\n",
+ argv[0], optarg, opt);
+ }
+ goto error;
+ }
+ }
+
+ if(!args->input) {
+ fprintf(stderr, "Missing input file.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ if(args->no_output) args->output = NULL;
+
+exit:
+ return res;
+error:
+ args_release(args);
+ goto exit;
+}
+
+/*******************************************************************************
+ * Grid header
+ ******************************************************************************/
+struct grid {
+ int32_t nx, ny, nz, ntimes; /* Definition */
+ int8_t is_z_irregular;
+ double lower[3]; /* Lower position of the grid */
+ double vxsz_x; /* Size of the voxel is X */
+ double vxsz_y; /* Size of the voxel in Y */
+ double* vxsz_z; /* Size of the voxel in Z */
+};
+
+#define GRID_NULL__ {0,0,0,0,0,{0,0,0},0,0,NULL}
+static const struct grid GRID_NULL = GRID_NULL__;
+
+static void
+grid_release(struct grid* grid)
+{
+ CHK(grid);
+ if(grid->vxsz_z) mem_rm(grid->vxsz_z);
+}
+
+/*******************************************************************************
+ * NetCDF helper functions and macros
+ ******************************************************************************/
+#define INVALID_ID -1
+
+#define NC(Func) { \
+ const int err__ = nc_ ## Func; \
+ if(err__ != NC_NOERR) { \
+ fprintf(stderr, "error:%i:%s\n", __LINE__, ncerr_to_str(err__)); \
+ abort(); \
+ } \
+ } (void)0
+
+static INLINE const char*
+ncerr_to_str(const int ncerr)
+{
+ const char* str = "NC_ERR_<UNKNOWN>";
+ switch(ncerr) {
+ case NC_EBADGRPID: str = "NC_EBADGRPID"; break;
+ case NC_EBADID: str = "NC_EBADID"; break;
+ case NC_EBADNAME: str = "NC_EBADNAME"; break;
+ case NC_ECHAR: str = "NC_ECHAR"; break;
+ case NC_EDIMMETA: str = "NC_EDIMMETA"; break;
+ case NC_EHDFERR: str = "NC_EHDFERR"; break;
+ case NC_ENOMEM: str = "NC_ENOMEM"; break;
+ case NC_ENOTATT: str = "NC_ENOTATT"; break;
+ case NC_ENOTVAR: str = "NC_ENOTVAR"; break;
+ case NC_ERANGE: str = "NC_ERANGE"; break;
+ case NC_NOERR: str = "NC_NOERR"; break;
+ }
+ return str;
+}
+
+static INLINE const char*
+nctype_to_str(const nc_type type)
+{
+ const char* str = "NC_TYPE_<UNKNOWN>";
+ switch(type) {
+ case NC_NAT: str = "NC_NAT"; break;
+ case NC_BYTE: str = "NC_BYTE"; break;
+ case NC_CHAR: str = "NC_CHAR"; break;
+ case NC_SHORT: str = "NC_SHORT"; break;
+ case NC_LONG: str = "NC_LONG"; break;
+ case NC_FLOAT: str = "NC_FLOAT"; break;
+ case NC_DOUBLE: str = "NC_DOUBLE"; break;
+ case NC_UBYTE: str = "NC_UBYTE"; break;
+ case NC_USHORT: str = "NC_USHORT"; break;
+ case NC_UINT: str = "NC_UINT"; break;
+ case NC_INT64: str = "NC_INT64"; break;
+ case NC_UINT64: str = "NC_UINT64"; break;
+ case NC_STRING: str = "NC_STRING"; break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+ return str;
+}
+
+static INLINE size_t
+sizeof_nctype(const nc_type type)
+{
+ size_t sz;
+ switch(type) {
+ case NC_BYTE:
+ case NC_CHAR:
+ case NC_UBYTE:
+ sz = 1; break;
+ case NC_SHORT:
+ case NC_USHORT:
+ sz = 2; break;
+ case NC_FLOAT:
+ case NC_INT:
+ case NC_UINT:
+ sz = 4; break;
+ case NC_DOUBLE:
+ case NC_INT64:
+ case NC_UINT64:
+ sz = 8; break;
+ default: FATAL("Unreachable cde.\n"); break;
+ }
+ return sz;
+}
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static res_T
+open_output_stream(const char* path, const int force_overwrite, FILE** stream)
+{
+ int fd = -1;
+ FILE* fp = NULL;
+ res_T res = RES_OK;
+ ASSERT(path);
+
+ if(force_overwrite) {
+ fp = fopen(path, "w");
+ if(!fp) {
+ fprintf(stderr, "Could not open the output file '%s'.\n", path);
+ goto error;
+ }
+ } else {
+ fd = open(path, O_CREAT|O_WRONLY|O_EXCL|O_TRUNC, S_IRUSR|S_IWUSR);
+ if(fd >= 0) {
+ fp = fdopen(fd, "w");
+ if(fp == NULL) {
+ fprintf(stderr, "Could not open the output file '%s'.\n", path);
+ goto error;
+ }
+ } else if(errno == EEXIST) {
+ fprintf(stderr,
+ "The output file '%s' already exists. Use -f to overwrite it.\n", path);
+ goto error;
+ } else {
+ fprintf(stderr,
+ "Unexpected error while opening the output file `'%s'.\n", path);
+ goto error;
+ }
+ }
+
+exit:
+ *stream = fp;
+ return res;
+error:
+ res = RES_IO_ERR;
+ if(fp) {
+ CHK(fclose(fp) == 0);
+ fp = NULL;
+ } else if(fd >= 0) {
+ CHK(close(fd) == 0);
+ }
+ goto exit;
+}
+
+
+static res_T
+setup_definition
+ (int nc, int32_t* nx, int32_t* ny, int32_t* nz, int32_t* ntimes)
+{
+ size_t len[4];
+ int dimids[4];
+ int id;
+ int ndims;
+ int err = NC_NOERR;
+ res_T res = RES_OK;
+ ASSERT(nx && ny && nz && ntimes);
+
+ err = nc_inq_varid(nc, "RCT", &id);
+ if(err != NC_NOERR) {
+ fprintf(stderr, "Could not inquire the 'RCT' variable -- %s\n",
+ ncerr_to_str(err));
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ NC(inq_varndims(nc, id, &ndims));
+ if(ndims != 4) {
+ fprintf(stderr, "The dimension of the 'RCT' variable must be 4.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ NC(inq_vardimid(nc, id, dimids));
+ NC(inq_dimlen(nc, dimids[0], len+0));
+ NC(inq_dimlen(nc, dimids[1], len+1));
+ NC(inq_dimlen(nc, dimids[2], len+2));
+ NC(inq_dimlen(nc, dimids[3], len+3));
+
+ *nx = (int32_t)len[3];
+ *ny = (int32_t)len[2];
+ *nz = (int32_t)len[1];
+ *ntimes = (int32_t)len[0];
+
+ if(*nx < 1 || *ny < 1 || *nz < 1 || *ntimes < 1) {
+ fprintf(stderr,
+ "The spatial and time definitions cannot be null.\n"
+ " #x = %i; #y = %i; #z = %i; #time = %i\n",
+ *nx, *ny, *nz, *ntimes);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
+static res_T
+setup_regular_dimension
+ (const int nc,
+ const char* var_name,
+ double* voxel_size,
+ double* lower_pos,
+ const int check)
+{
+ char* mem = NULL;
+ size_t len, start, i;
+ double vxsz;
+ nc_type type;
+ int varid, dimid, ndims;
+ int err = NC_NOERR;
+ res_T res = RES_OK;
+ ASSERT(nc && var_name && voxel_size && lower_pos);
+
+ err = nc_inq_varid(nc, var_name, &varid);
+ if(err != NC_NOERR) {
+ fprintf(stderr, "Could not inquire the '%s' variable -- %s\n",
+ var_name, ncerr_to_str(err));
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ NC(inq_varndims(nc, varid, &ndims));
+ if(ndims != 1) {
+ fprintf(stderr, "The dimension of the '%s' variable must be 1 .\n", var_name);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ NC(inq_vardimid(nc, varid, &dimid));
+ NC(inq_dimlen(nc, dimid, &len));
+ ASSERT(len > 0);
+
+ NC(inq_vartype(nc, varid, &type));
+ if(type != NC_DOUBLE && type != NC_FLOAT) {
+ fprintf(stderr, "The type of the '%s' variable must be float or double.\n",
+ var_name);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ #define AT(Mem, Id) (type == NC_FLOAT ? ((float*)Mem)[Id] : ((double*)Mem)[Id])
+ if(!check) {
+ char* raw = alloca(8);
+ start = 0;
+ len = 1;
+ NC(get_vara(nc, varid, &start, &len, raw));
+ vxsz = AT(raw, 0) * 2; /* Assume that the grid starts at 0 */
+ } else {
+ mem = mem_alloc(len*sizeof_nctype(type));
+ if(!mem) {
+ fprintf(stderr, "Could not allocate memory for the '%s' variable.\n", var_name);
+ res = RES_MEM_ERR;
+ goto error;
+ }
+ start = 0;
+ NC(get_vara(nc, varid, &start, &len, mem));
+ vxsz = AT(mem, 0) * 2; /* Assume that the grid starts from 0 */
+ }
+
+ if(vxsz <= 0) {
+ fprintf(stderr, "The '%s' variable can't have negative or null values.\n",
+ var_name);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ if(check) {
+ /* Check that the dimension is regular */
+ FOR_EACH(i, 1, len) {
+ if(!eq_eps(AT(mem, i) - AT(mem, i-1), vxsz, 1.e-6)) {
+ fprintf(stderr, "The voxel size of the '%s' variable must be regular.\n",
+ var_name);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ }
+ }
+ #undef AT
+
+ *voxel_size = vxsz;
+ *lower_pos = 0;
+exit:
+ if(mem) mem_rm(mem);
+ return res;
+error:
+ goto exit;
+}
+
+/* Return 1 if Z is irregular */
+static int
+setup_Z_dimension
+ (const int nc,
+ double** voxel_size,
+ double* lower_pos,
+ int8_t* is_z_irregular,
+ const int check)
+{
+ enum { Z, X, Y, NDIMS }; /* Helper constants */
+ char* mem = NULL; /* Memory where NetCDF data are copied */
+ char* mem2 = NULL; /* Temporary memory used to check the VLEV data */
+ size_t len[NDIMS];
+ size_t start[NDIMS];
+ size_t dim_len[NDIMS];
+ size_t z; /* Id over the Z dimension */
+ size_t i; /* Common id */
+ double vxsz; /* Voxel size */
+ double* pvxsz = NULL; /* Pointer toward voxel sizes */
+ nc_type type; /* Type of the NetCDF VLEV variable */
+ int varid; /* NetCDF id of the VLEV variable */
+ int dimids[NDIMS]; /* NetCDF id of the VLEV dimensions */
+ int irregular; /* Boolean defining if the Z dimension is irregular */
+ int ndims; /* #dims of the VLEV variable */
+ int err = NC_NOERR;
+ res_T res = RES_OK;
+ ASSERT(nc && voxel_size && lower_pos && is_z_irregular);
+
+ /* Retrieve the VLEV variable */
+ err = nc_inq_varid(nc, "VLEV", &varid);
+ if(err != NC_NOERR) {
+ fprintf(stderr, "Could not inquire the 'VLEV' variable -- %s\n",
+ ncerr_to_str(err));
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Inquire the VLEV dimensionality */
+ NC(inq_varndims(nc, varid, &ndims));
+ if(ndims != NDIMS) {
+ fprintf(stderr, "The dimension of the 'VLEV' variable must be %i .\n",
+ NDIMS);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ NC(inq_vardimid(nc, varid, dimids));
+ FOR_EACH(i, 0, NDIMS) NC(inq_dimlen(nc, dimids[i], len+i));
+
+ /* Inquire the type of the VLEV variable */
+ NC(inq_vartype(nc, varid, &type));
+ if(type != NC_DOUBLE && type != NC_FLOAT) {
+ fprintf(stderr, "The type of 'VLEV' variable must be float or double.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Allocate te memory space where a Z column of the VLEV variable is going to
+ * be copied */
+ mem = mem_alloc(len[Z]*sizeof_nctype(type));
+ if(!mem) {
+ fprintf(stderr, "Could not allocate memory for the 'VLEV' variable.\n");
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
+ /* Get only one Z column */
+ start[X] = 0; dim_len[X] = 1;
+ start[Y] = 0; dim_len[Y] = 1;
+ start[Z] = 0; dim_len[Z] = len[Z];
+ NC(get_vara(nc, varid, start, dim_len, mem));
+
+ #define AT(Mem, Id) (type == NC_FLOAT ? ((float*)Mem)[Id] : ((double*)Mem)[Id])
+
+ /* Assume that the Z dimension starts from 0 */
+ vxsz = AT(mem, 0) * 2.0;
+ if(vxsz <= 0) {
+ fprintf(stderr, "The 'VLEV' variable can't have negative or null values.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Check if the Z dimension is regular */
+ FOR_EACH(z, 1, len[Z]) {
+ if(!eq_eps(AT(mem, z) - AT(mem, z-1), vxsz, 1.e-6)) break;
+ }
+ irregular = (z != len[Z]);
+
+ /* Setup the size of the voxel */
+ pvxsz = mem_alloc(sizeof(double) * (irregular ? len[Z] : 1));
+ if(!pvxsz) {
+ fprintf(stderr,
+ "Could not allocate memory for the voxel size along the Z dimension.\n");
+ res = RES_MEM_ERR;
+ goto error;
+ }
+ pvxsz[0] = vxsz;
+
+ if(irregular) {
+ FOR_EACH(z, 1, len[Z]) {
+ pvxsz[z] = (AT(mem, z) - (AT(mem, z-1) + pvxsz[z-1]*0.5)) * 2.0;
+ if(pvxsz[z] <= 0) {
+ fprintf(stderr,
+ "The 'VLEV' variable can't have negative or null values.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ }
+ }
+
+ /* Check that the others columns have the same voxel size */
+ if(check) {
+ size_t ix, iy;
+
+ mem2 = mem_alloc(len[Z]*sizeof_nctype(type));
+ if(!mem2) {
+ fprintf(stderr,
+ "Could not allocate memory to check the 'VLEV' variable.\n");
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
+ FOR_EACH(iy, 0, len[Y]) {
+ start[Y] = iy;
+ FOR_EACH(ix, 1, len[X]) {
+ start[X] = ix;
+ NC(get_vara(nc, varid, start, dim_len, mem2));
+ FOR_EACH(z, 0, len[Z]) {
+ if(!eq_eps(AT(mem2, z), AT(mem, z), 1.e-6)) {
+ fprintf(stderr,
+ "The Z columns of the 'VLEV' variable must be the same.\n");
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ }
+ }
+ }
+ }
+ #undef AT
+
+ *lower_pos = 0;
+ *voxel_size = pvxsz;
+ *is_z_irregular = (int8_t)irregular;
+
+exit:
+ if(mem2) mem_rm(mem2);
+ if(mem) mem_rm(mem);
+ return res;
+error:
+ goto exit;
+}
+
+static res_T
+write_data(int nc, const char* var_name, FILE* stream)
+{
+ enum { TIME, Z, X, Y, NDIMS }; /* Helper constants */
+
+ char* mem = NULL; /* Memory where NetCDF data are copied */
+ char* mem2 = NULL; /* Tmp memory where NetCDF data are converted */
+ size_t start[NDIMS];
+ size_t dim_len[NDIMS];
+ size_t len[NDIMS];
+ size_t slice_len; /* #elements per slice the 3D grid */
+ size_t i; /* Common id */
+ size_t t; /* Id over the time dimension */
+ size_t z; /* Id over the Z dimension */
+ nc_type type;
+ int varid; /* NetCDF id of the data to write */
+ int ndims; /* #dimensions of the NetCDF data to write */
+ int dimids[NDIMS]; /* NetCDF id of the data dimension */
+ int err = NC_NOERR;
+ res_T res = RES_OK;
+ CHK(var_name);
+
+ /* Retrieve the NetCDF id of the variable */
+ err = nc_inq_varid(nc, var_name, &varid);
+ if(err != NC_NOERR) {
+ fprintf(stderr, "Could not inquire the '%s' variable -- %s\n",
+ var_name, ncerr_to_str(err));
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Inquire the dimensions of the variable */
+ NC(inq_varndims(nc, varid, &ndims));
+ if(ndims != NDIMS) {
+ fprintf(stderr, "The dimension of the '%s' variable must be %i .\n",
+ var_name, NDIMS);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ NC(inq_vardimid(nc, varid, dimids));
+ FOR_EACH(i, 0, NDIMS) NC(inq_dimlen(nc, dimids[i], len+i));
+
+ /* Inquire the type of the variable */
+ NC(inq_vartype(nc, varid, &type));
+ if(type != NC_DOUBLE && type != NC_FLOAT) {
+ fprintf(stderr,
+ "The type of the '%s' variable cannot be %s. Expecting floating point data.\n",
+ var_name, nctype_to_str(type));
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Alloc the memory space where the variable data wille be copied. Note that,
+ * in order to reduce memory consumption, only one slice is read at a given
+ * time and Z position. */
+ slice_len = len[X]*len[Y]; /* # elements per slice */
+ mem = mem_alloc(slice_len*sizeof_nctype(type));
+
+ /* Allocate a temporary memory space to store the variable data converted in
+ * double precision if they are stored in single precision. */
+ mem2 = type==NC_DOUBLE ? mem : mem_alloc(slice_len*sizeof_nctype(NC_DOUBLE));
+
+ if(!mem || !mem2) {
+ fprintf(stderr, "Could not allocate memory for the '%s' variable.\n",
+ var_name);
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
+ /* Read the whole slice */
+ start[X] = 0; dim_len[X] = len[X];
+ start[Y] = 0; dim_len[Y] = len[Y];
+
+ FOR_EACH(t, 0, len[TIME]) {
+ start[TIME] = t, dim_len[TIME] = 1;
+ FOR_EACH(z, 0, len[Z]) {
+ size_t n;
+ start[Z] = z; dim_len[Z] = 1;
+
+ NC(get_vara(nc, varid, start, dim_len, mem)); /* Fetch the slice data */
+
+ /* In order to be compliant with the htgop fileformat, single precision
+ * data must be converted in double precision */
+ if(type == NC_FLOAT) {
+ FOR_EACH(i, 0, slice_len) ((double*)mem2)[i] = ((float*)mem)[i];
+ }
+
+ /* Write data */
+ n = fwrite(mem2, sizeof_nctype(NC_DOUBLE), slice_len, stream);
+ if(n != slice_len) {
+ fprintf(stderr, "Error writing data of the '%s' variable -- '%s'.\n",
+ var_name, strerror(errno));
+ res = RES_IO_ERR;
+ goto error;
+ }
+ }
+ }
+
+exit:
+ if(mem2 != mem) mem_rm(mem2);
+ if(mem) mem_rm(mem);
+ return res;
+error:
+ goto exit;
+}
+
+/*******************************************************************************
+ * Program
+ ******************************************************************************/
+int
+main(int argc, char** argv)
+{
+ FILE* stream = stdout;
+ struct grid grid = GRID_NULL;
+ struct args args = ARGS_DEFAULT;
+ int err = 0;
+ int err_nc = NC_NOERR;
+ int nc = INVALID_ID;
+ res_T res = RES_OK;
+
+ res = args_init(&args, argc, argv);
+ if(res != RES_OK) goto error;
+ if(args.quit) goto exit;
+
+ if(args.output) {
+ res = open_output_stream(args.output, args.force_overwrite, &stream);
+ if(res != RES_OK) goto error;
+ }
+
+ err_nc = nc_open(args.input, NC_WRITE|NC_MMAP, &nc);
+ if(err_nc != NC_NOERR) {
+ fprintf(stderr, "Error opening file `%s' -- %s.\n",
+ args.input, ncerr_to_str(err_nc));
+ goto exit;
+ }
+
+ #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0
+ CALL(setup_definition(nc, &grid.nx, &grid.ny, &grid.nz, &grid.ntimes));
+ CALL(setup_regular_dimension
+ (nc, "W_E_direction", &grid.vxsz_x, grid.lower+0, args.check));
+ CALL(setup_regular_dimension
+ (nc, "S_N_direction", &grid.vxsz_y, grid.lower+1, args.check));
+ CALL(setup_Z_dimension
+ (nc, &grid.vxsz_z, grid.lower+2, &grid.is_z_irregular, args.check));
+ #undef CALL
+
+ #define WRITE(Var, N, Name) { \
+ if(fwrite((Var), sizeof(*(Var)), (N), stream) != (N)) { \
+ fprintf(stderr, "Error writing the %s.\n", (Name)); \
+ res = RES_IO_ERR; \
+ goto error; \
+ } \
+ } (void)0
+ if(!args.no_output) {
+ const int64_t pagesize = (int64_t)args.pagesize;
+
+ WRITE(&pagesize, 1, "page size");
+ WRITE(&grid.is_z_irregular, 1, "'irregular' Z boolean");
+ WRITE(&grid.nx, 1, "X definition");
+ WRITE(&grid.ny, 1, "Y definition");
+ WRITE(&grid.nz, 1, "Z definition");
+ WRITE(&grid.ntimes, 1, "time definition");
+ WRITE(grid.lower, 3, "lower position");
+ WRITE(&grid.vxsz_x, 1, "X voxel size");
+ WRITE(&grid.vxsz_y, 1, "Y voxel size");
+ WRITE(grid.vxsz_z, grid.is_z_irregular?(size_t)grid.nz:1, "Z voxel size(s)");
+ /* padding */
+ fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
+ write_data(nc, "RVT", stream);
+ /* padding */
+ fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
+ write_data(nc, "RCT", stream);
+ /*padding */
+ fseek(stream, ALIGN_SIZE(ftell(stream), (off_t)pagesize), SEEK_SET);
+ }
+ #undef WRITE
+
+exit:
+ grid_release(&grid);
+ if(nc != INVALID_ID) NC(close(nc));
+ if(mem_allocated_size() != 0) {
+ fprintf(stderr, "Memory leaks: %lu Bytes\n",
+ (unsigned long)mem_allocated_size());
+ err = -1;
+ }
+ return err;
+error:
+ err = -1;
+ goto exit;
+}
+
diff --git a/src/les2htcop.h.in b/src/les2htcp.h.in
diff --git a/src/test_htcop.c b/src/test_htcop.c
@@ -1,65 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#include "htcop.h"
-#include "test_htcop_utils.h"
-
-#include <rsys/logger.h>
-
-static void
-log_stream(const char* msg, void* ctx)
-{
- ASSERT(msg);
- (void)msg, (void)ctx;
- printf("%s\n", msg);
-}
-
-int
-main(int argc, char** argv)
-{
- struct mem_allocator allocator;
- struct logger logger;
- struct htcop* htcop = NULL;
- (void)argc, (void)argv;
-
- CHK(htcop_create(NULL, NULL, 0, NULL) == RES_BAD_ARG);
- CHK(htcop_create(NULL, NULL, 0, &htcop) == RES_OK);
-
- CHK(htcop_ref_get(NULL) == RES_BAD_ARG);
- CHK(htcop_ref_get(htcop) == RES_OK);
- CHK(htcop_ref_put(NULL) == RES_BAD_ARG);
- CHK(htcop_ref_put(htcop) == RES_OK);
- CHK(htcop_ref_put(htcop) == RES_OK);
-
- CHK(mem_init_proxy_allocator(&allocator, &mem_default_allocator) == RES_OK);
- CHK(htcop_create(NULL, &allocator, 1, &htcop) == RES_OK);
- CHK(htcop_ref_put(htcop) == RES_OK);
-
- CHK(logger_init(&allocator, &logger) == RES_OK);
- logger_set_stream(&logger, LOG_OUTPUT, log_stream, NULL);
- logger_set_stream(&logger, LOG_ERROR, log_stream, NULL);
- logger_set_stream(&logger, LOG_WARNING, log_stream, NULL);
-
- CHK(htcop_create(&logger, &allocator, 0, &htcop) == RES_OK);
- CHK(htcop_ref_put(htcop) == RES_OK);
- CHK(htcop_create(&logger, NULL, 0, &htcop) == RES_OK);
- CHK(htcop_ref_put(htcop) == RES_OK);
-
- logger_release(&logger);
- check_memory_allocator(&allocator);
- mem_shutdown_proxy_allocator(&allocator);
- CHK(mem_allocated_size() == 0);
- return 0;
-}
diff --git a/src/test_htcop_load.c b/src/test_htcop_load.c
@@ -1,120 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#include "htcop.h"
-#include "test_htcop_utils.h"
-
-#include <unistd.h>
-
-int
-main(int argc, char** argv)
-{
- int64_t pagesize;
- struct htcop* htcop = NULL;
- struct htcop_desc desc = HTCOP_DESC_NULL;
- FILE* stream = NULL;
- int8_t i8;
- int32_t i32[4];
- double dbl[3];
- size_t i, n;
- size_t x, y, z, t;
- (void)argc, (void)argv;
-
- stream = tmpfile();
-
- pagesize = (int64_t)sysconf(_SC_PAGESIZE); /* pagesize */
- CHK(fwrite(&pagesize, sizeof(pagesize), 1, stream) == 1);
- i8 = 0; /* irregular Z */
- CHK(fwrite(&i8, sizeof(i8), 1, stream) == 1);
- i32[0] = 2; /* Definition */
- i32[1] = 2;
- i32[2] = 3;
- i32[3] = 1;
- CHK(fwrite(i32, sizeof(i32[0]), 4, stream) == 4);
- dbl[0] = 0; /* Lower position */
- dbl[1] = 0;
- dbl[2] = 0;
- CHK(fwrite(dbl, sizeof(dbl[0]), 3, stream) == 3);
- dbl[0] = 1; /* Voxel size */
- dbl[1] = 2;
- dbl[2] = 3;
- CHK(fwrite(dbl, sizeof(dbl[0]), 3, stream) == 3);
-
- fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
- FOR_EACH(i, 0, i32[0]*i32[1]*i32[2]*i32[3]) { /* RVT */
- dbl[0] = (double)i;
- CHK(fwrite(dbl, sizeof(dbl[0]), 1, stream) == 1);
- }
-
- fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
- FOR_EACH(i, 0, i32[0]*i32[1]*i32[2]*i32[3]) { /* RCT */
- dbl[0] =-(double)i;
- CHK(fwrite(dbl, sizeof(dbl[0]), 1, stream) == 1);
- }
- fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
-
- rewind(stream);
-
- CHK(htcop_create(NULL, &mem_default_allocator, 1, &htcop) == RES_OK);
-
- CHK(htcop_get_desc(NULL, &desc) == RES_BAD_ARG);
- CHK(htcop_get_desc(htcop, NULL) == RES_BAD_ARG);
- CHK(htcop_get_desc(htcop, &desc) == RES_BAD_ARG);
-
- CHK(htcop_load_stream(NULL, stream) == RES_BAD_ARG);
- CHK(htcop_load_stream(htcop, NULL) == RES_BAD_ARG);
- CHK(htcop_load_stream(htcop, stream) == RES_OK);
- fclose(stream);
-
- CHK(htcop_get_desc(NULL, &desc) == RES_BAD_ARG);
- CHK(htcop_get_desc(htcop, NULL) == RES_BAD_ARG);
- CHK(htcop_get_desc(htcop, &desc) == RES_OK);
-
- CHK(desc.pagesize == 4096);
- CHK(desc.irregular_z == 0);
- CHK(desc.spatial_definition[0] == 2);
- CHK(desc.spatial_definition[1] == 2);
- CHK(desc.spatial_definition[2] == 3);
- CHK(desc.time_definition == 1);
- CHK(desc.lower[0] == 0);
- CHK(desc.lower[1] == 0);
- CHK(desc.lower[2] == 0);
- CHK(desc.vxsz_x == 1);
- CHK(desc.vxsz_y == 2);
- CHK(desc.vxsz_z[0] == 3);
-
- n = desc.spatial_definition[0]
- * desc.spatial_definition[1]
- * desc.spatial_definition[2]
- * desc.time_definition;
- FOR_EACH(i, 0,n) {
- CHK(desc.RVT[i] == (double)i);
- CHK(desc.RCT[i] ==-(double)i);
- }
- i = 0;
- FOR_EACH(t, 0, desc.time_definition) {
- FOR_EACH(z, 0, desc.spatial_definition[2]) {
- FOR_EACH(y, 0, desc.spatial_definition[1]) {
- FOR_EACH(x, 0, desc.spatial_definition[0]) {
- CHK(htcop_desc_RVT_at(&desc, x, y, z, t) == (double)i);
- CHK(htcop_desc_RCT_at(&desc, x, y, z, t) ==-(double)i);
- ++i;
- }}}}
- CHK(htcop_ref_put(htcop) == RES_OK);
-
- CHK(mem_allocated_size() == 0);
- return 0;
-}
-
diff --git a/src/test_htcop_load_from_file.c b/src/test_htcop_load_from_file.c
@@ -1,75 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#include "htcop.h"
-#include "test_htcop_utils.h"
-
-int
-main(int argc, char** argv)
-{
- struct htcop* htcop = NULL;
- struct htcop_desc desc = HTCOP_DESC_NULL;
- size_t i, n;
-
- if(argc < 2) {
- fprintf(stderr, "Usage: %s HTCOP_FILE\n", argv[0]);
- return -1;
- }
-
- CHK(htcop_create(NULL, &mem_default_allocator, 1, &htcop) == RES_OK);
- CHK(htcop_load(htcop, argv[1]) == RES_OK);
-
- CHK(htcop_get_desc(htcop, &desc) == RES_OK);
-
- printf("pagesize: %lu\n", (unsigned long)desc.pagesize);
- printf("irregular Z: %i\n", desc.irregular_z);
- printf("#X: %lu; #Y: %lu; #Z: %lu; #times: %lu\n",
- desc.spatial_definition[0],
- desc.spatial_definition[1],
- desc.spatial_definition[2],
- desc.time_definition);
- printf("lower pos: %g %g %g\n", SPLIT3(desc.lower));
- printf("voxel size: %g %g ", desc.vxsz_x, desc.vxsz_y);
- if(!desc.irregular_z) {
- printf("%g\n", desc.vxsz_z[0]);
- } else {
- printf("{");
- FOR_EACH(i, 0, desc.spatial_definition[2]) {
- printf("%g", desc.vxsz_z[i]);
- if(i != desc.spatial_definition[2]-1) printf(", ");
- }
- printf("}");
- }
-
- n = desc.spatial_definition[0]
- * desc.spatial_definition[1]
- * desc.spatial_definition[2]
- * desc.time_definition;
-
- printf("RVT:\n");
- FOR_EACH(i, 0, n) {
- printf("%g\n", desc.RVT[i]);
- }
-
- printf("RCT:\n");
- FOR_EACH(i, 0, n) {
- printf("%g\n", desc.RCT[i]);
- }
-
- CHK(htcop_ref_put(htcop) == RES_OK);
- check_memory_allocator(&mem_default_allocator);
- CHK(mem_allocated_size() == 0);
- return 0;
-}
diff --git a/src/test_htcop_utils.h b/src/test_htcop_utils.h
@@ -1,34 +0,0 @@
-/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>. */
-
-#ifndef TEST_HTCOP_UTILS_H
-#define TEST_HTCOP_UTILS_H
-
-#include <rsys/mem_allocator.h>
-#include <stdio.h>
-
-static INLINE void
-check_memory_allocator(struct mem_allocator* allocator)
-{
- if(MEM_ALLOCATED_SIZE(allocator)) {
- char dump[512];
- MEM_DUMP(allocator, dump, sizeof(dump)/sizeof(char));
- fprintf(stderr, "%s\n", dump);
- FATAL("Memory leaks\n");
- }
-}
-
-#endif /* TEST_HTCOP_UTILS_H */
-
diff --git a/src/test_htcp.c b/src/test_htcp.c
@@ -0,0 +1,65 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "htcp.h"
+#include "test_htcp_utils.h"
+
+#include <rsys/logger.h>
+
+static void
+log_stream(const char* msg, void* ctx)
+{
+ ASSERT(msg);
+ (void)msg, (void)ctx;
+ printf("%s\n", msg);
+}
+
+int
+main(int argc, char** argv)
+{
+ struct mem_allocator allocator;
+ struct logger logger;
+ struct htcp* htcp = NULL;
+ (void)argc, (void)argv;
+
+ CHK(htcp_create(NULL, NULL, 0, NULL) == RES_BAD_ARG);
+ CHK(htcp_create(NULL, NULL, 0, &htcp) == RES_OK);
+
+ CHK(htcp_ref_get(NULL) == RES_BAD_ARG);
+ CHK(htcp_ref_get(htcp) == RES_OK);
+ CHK(htcp_ref_put(NULL) == RES_BAD_ARG);
+ CHK(htcp_ref_put(htcp) == RES_OK);
+ CHK(htcp_ref_put(htcp) == RES_OK);
+
+ CHK(mem_init_proxy_allocator(&allocator, &mem_default_allocator) == RES_OK);
+ CHK(htcp_create(NULL, &allocator, 1, &htcp) == RES_OK);
+ CHK(htcp_ref_put(htcp) == RES_OK);
+
+ CHK(logger_init(&allocator, &logger) == RES_OK);
+ logger_set_stream(&logger, LOG_OUTPUT, log_stream, NULL);
+ logger_set_stream(&logger, LOG_ERROR, log_stream, NULL);
+ logger_set_stream(&logger, LOG_WARNING, log_stream, NULL);
+
+ CHK(htcp_create(&logger, &allocator, 0, &htcp) == RES_OK);
+ CHK(htcp_ref_put(htcp) == RES_OK);
+ CHK(htcp_create(&logger, NULL, 0, &htcp) == RES_OK);
+ CHK(htcp_ref_put(htcp) == RES_OK);
+
+ logger_release(&logger);
+ check_memory_allocator(&allocator);
+ mem_shutdown_proxy_allocator(&allocator);
+ CHK(mem_allocated_size() == 0);
+ return 0;
+}
diff --git a/src/test_htcp_load.c b/src/test_htcp_load.c
@@ -0,0 +1,120 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "htcp.h"
+#include "test_htcp_utils.h"
+
+#include <unistd.h>
+
+int
+main(int argc, char** argv)
+{
+ int64_t pagesize;
+ struct htcp* htcp = NULL;
+ struct htcp_desc desc = HTCP_DESC_NULL;
+ FILE* stream = NULL;
+ int8_t i8;
+ int32_t i32[4];
+ double dbl[3];
+ size_t i, n;
+ size_t x, y, z, t;
+ (void)argc, (void)argv;
+
+ stream = tmpfile();
+
+ pagesize = (int64_t)sysconf(_SC_PAGESIZE); /* pagesize */
+ CHK(fwrite(&pagesize, sizeof(pagesize), 1, stream) == 1);
+ i8 = 0; /* irregular Z */
+ CHK(fwrite(&i8, sizeof(i8), 1, stream) == 1);
+ i32[0] = 2; /* Definition */
+ i32[1] = 2;
+ i32[2] = 3;
+ i32[3] = 1;
+ CHK(fwrite(i32, sizeof(i32[0]), 4, stream) == 4);
+ dbl[0] = 0; /* Lower position */
+ dbl[1] = 0;
+ dbl[2] = 0;
+ CHK(fwrite(dbl, sizeof(dbl[0]), 3, stream) == 3);
+ dbl[0] = 1; /* Voxel size */
+ dbl[1] = 2;
+ dbl[2] = 3;
+ CHK(fwrite(dbl, sizeof(dbl[0]), 3, stream) == 3);
+
+ fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
+ FOR_EACH(i, 0, i32[0]*i32[1]*i32[2]*i32[3]) { /* RVT */
+ dbl[0] = (double)i;
+ CHK(fwrite(dbl, sizeof(dbl[0]), 1, stream) == 1);
+ }
+
+ fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
+ FOR_EACH(i, 0, i32[0]*i32[1]*i32[2]*i32[3]) { /* RCT */
+ dbl[0] =-(double)i;
+ CHK(fwrite(dbl, sizeof(dbl[0]), 1, stream) == 1);
+ }
+ fseek(stream, ALIGN_SIZE(ftell(stream), pagesize), SEEK_SET); /* padding */
+
+ rewind(stream);
+
+ CHK(htcp_create(NULL, &mem_default_allocator, 1, &htcp) == RES_OK);
+
+ CHK(htcp_get_desc(NULL, &desc) == RES_BAD_ARG);
+ CHK(htcp_get_desc(htcp, NULL) == RES_BAD_ARG);
+ CHK(htcp_get_desc(htcp, &desc) == RES_BAD_ARG);
+
+ CHK(htcp_load_stream(NULL, stream) == RES_BAD_ARG);
+ CHK(htcp_load_stream(htcp, NULL) == RES_BAD_ARG);
+ CHK(htcp_load_stream(htcp, stream) == RES_OK);
+ fclose(stream);
+
+ CHK(htcp_get_desc(NULL, &desc) == RES_BAD_ARG);
+ CHK(htcp_get_desc(htcp, NULL) == RES_BAD_ARG);
+ CHK(htcp_get_desc(htcp, &desc) == RES_OK);
+
+ CHK(desc.pagesize == 4096);
+ CHK(desc.irregular_z == 0);
+ CHK(desc.spatial_definition[0] == 2);
+ CHK(desc.spatial_definition[1] == 2);
+ CHK(desc.spatial_definition[2] == 3);
+ CHK(desc.time_definition == 1);
+ CHK(desc.lower[0] == 0);
+ CHK(desc.lower[1] == 0);
+ CHK(desc.lower[2] == 0);
+ CHK(desc.vxsz_x == 1);
+ CHK(desc.vxsz_y == 2);
+ CHK(desc.vxsz_z[0] == 3);
+
+ n = desc.spatial_definition[0]
+ * desc.spatial_definition[1]
+ * desc.spatial_definition[2]
+ * desc.time_definition;
+ FOR_EACH(i, 0,n) {
+ CHK(desc.RVT[i] == (double)i);
+ CHK(desc.RCT[i] ==-(double)i);
+ }
+ i = 0;
+ FOR_EACH(t, 0, desc.time_definition) {
+ FOR_EACH(z, 0, desc.spatial_definition[2]) {
+ FOR_EACH(y, 0, desc.spatial_definition[1]) {
+ FOR_EACH(x, 0, desc.spatial_definition[0]) {
+ CHK(htcp_desc_RVT_at(&desc, x, y, z, t) == (double)i);
+ CHK(htcp_desc_RCT_at(&desc, x, y, z, t) ==-(double)i);
+ ++i;
+ }}}}
+ CHK(htcp_ref_put(htcp) == RES_OK);
+
+ CHK(mem_allocated_size() == 0);
+ return 0;
+}
+
diff --git a/src/test_htcp_load_from_file.c b/src/test_htcp_load_from_file.c
@@ -0,0 +1,75 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "htcp.h"
+#include "test_htcp_utils.h"
+
+int
+main(int argc, char** argv)
+{
+ struct htcp* htcp = NULL;
+ struct htcp_desc desc = HTCP_DESC_NULL;
+ size_t i, n;
+
+ if(argc < 2) {
+ fprintf(stderr, "Usage: %s HTCP_FILE\n", argv[0]);
+ return -1;
+ }
+
+ CHK(htcp_create(NULL, &mem_default_allocator, 1, &htcp) == RES_OK);
+ CHK(htcp_load(htcp, argv[1]) == RES_OK);
+
+ CHK(htcp_get_desc(htcp, &desc) == RES_OK);
+
+ printf("pagesize: %lu\n", (unsigned long)desc.pagesize);
+ printf("irregular Z: %i\n", desc.irregular_z);
+ printf("#X: %lu; #Y: %lu; #Z: %lu; #times: %lu\n",
+ desc.spatial_definition[0],
+ desc.spatial_definition[1],
+ desc.spatial_definition[2],
+ desc.time_definition);
+ printf("lower pos: %g %g %g\n", SPLIT3(desc.lower));
+ printf("voxel size: %g %g ", desc.vxsz_x, desc.vxsz_y);
+ if(!desc.irregular_z) {
+ printf("%g\n", desc.vxsz_z[0]);
+ } else {
+ printf("{");
+ FOR_EACH(i, 0, desc.spatial_definition[2]) {
+ printf("%g", desc.vxsz_z[i]);
+ if(i != desc.spatial_definition[2]-1) printf(", ");
+ }
+ printf("}");
+ }
+
+ n = desc.spatial_definition[0]
+ * desc.spatial_definition[1]
+ * desc.spatial_definition[2]
+ * desc.time_definition;
+
+ printf("RVT:\n");
+ FOR_EACH(i, 0, n) {
+ printf("%g\n", desc.RVT[i]);
+ }
+
+ printf("RCT:\n");
+ FOR_EACH(i, 0, n) {
+ printf("%g\n", desc.RCT[i]);
+ }
+
+ CHK(htcp_ref_put(htcp) == RES_OK);
+ check_memory_allocator(&mem_default_allocator);
+ CHK(mem_allocated_size() == 0);
+ return 0;
+}
diff --git a/src/test_htcp_utils.h b/src/test_htcp_utils.h
@@ -0,0 +1,34 @@
+/* Copyright (C) 2018 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#ifndef TEST_HTCP_UTILS_H
+#define TEST_HTCP_UTILS_H
+
+#include <rsys/mem_allocator.h>
+#include <stdio.h>
+
+static INLINE void
+check_memory_allocator(struct mem_allocator* allocator)
+{
+ if(MEM_ALLOCATED_SIZE(allocator)) {
+ char dump[512];
+ MEM_DUMP(allocator, dump, sizeof(dump)/sizeof(char));
+ fprintf(stderr, "%s\n", dump);
+ FATAL("Memory leaks\n");
+ }
+}
+
+#endif /* TEST_HTCP_UTILS_H */
+
