diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..e371feb
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,3 @@
+perf.data*
+build/
+.vscode/
diff --git a/CMakeLists.txt b/CMakeLists.txt
new file mode 100644
index 0000000..c1b16a2
--- /dev/null
+++ b/CMakeLists.txt
@@ -0,0 +1,63 @@
+cmake_minimum_required (VERSION 3.9)
+cmake_policy(SET CMP0077 NEW)
+
+project(meshed
+ VERSION 1.0.0
+ LANGUAGES C
+)
+
+option(USE_ASAN "Enable or disable address sanitizer" OFF)
+option(USE_USAN "Enable or disable undefined sanitizer" OFF)
+
+find_package(GDAL REQUIRED)
+find_package(OpenMP REQUIRED)
+
+set(CMAKE_C_STANDARD 99)
+
+set(WARNING_FLAGS
+ -Waddress
+ -Wall
+ -Wdeprecated
+ -Wextra
+ -Wimplicit-fallthrough
+ -Wmisleading-indentation
+ -Wmissing-prototypes
+ -Wshadow
+ -Wshift-sign-overflow
+ -Wswitch
+ -Wuninitialized
+ -Wunreachable-code
+ -Wunused
+ -Wunused-but-set-variable
+ -Wunused-function
+ # -Werror
+)
+
+add_executable(meshed.exe
+ delineate.c
+ delineate_lessmem.c
+ delineate_moremem.c
+ gettimeofday.c
+ hierarchy.c
+ main.c
+ outlet_list.c
+ outlets.c
+ raster.c
+ timeval_diff.c
+)
+
+target_link_libraries(meshed.exe PRIVATE GDAL::GDAL OpenMP::OpenMP_C m)
+target_compile_options(meshed.exe PRIVATE ${WARNING_FLAGS})
+set_property(TARGET meshed.exe PROPERTY INTERPROCEDURAL_OPTIMIZATION TRUE)
+if (CMAKE_C_COMPILER_ID STREQUAL "GNU")
+ target_compile_options(meshed.exe PRIVATE -foptimize-sibling-calls)
+endif()
+if(USE_ASAN)
+ target_compile_options(meshed.exe PRIVATE -fsanitize=address)
+ target_link_options(meshed.exe PRIVATE -fsanitize=address)
+endif()
+
+if(USE_USAN)
+ target_compile_options(meshed.exe PRIVATE -fsanitize=undefined)
+ target_link_options(meshed.exe PRIVATE -fsanitize=undefined)
+endif()
diff --git a/Makefile b/Makefile
deleted file mode 100644
index 4ebdff7..0000000
--- a/Makefile
+++ /dev/null
@@ -1,29 +0,0 @@
-ifeq ($(OS),Windows_NT)
- GDAL_CFLAGS=-I/c/OSGeo4W/include
- GDAL_LIBS=/c/OSGeo4W/lib/gdal_i.lib
- EXT=.exe
-else
- GDAL_LIBS=`gdal-config --libs`
-endif
-CFLAGS=-Wall -Werror -O3 -fopenmp $(GDAL_CFLAGS)
-LDFLAGS=-O3 -fopenmp -lm
-
-all: meshed$(EXT)
-
-clean:
- $(RM) *.o
-
-meshed$(EXT): \
- main.o \
- timeval_diff.o \
- raster.o \
- outlet_list.o \
- outlets.o \
- delineate.o \
- delineate_lessmem.o \
- delineate_moremem.o \
- hierarchy.o
- $(CC) $(LDFLAGS) -o $@ $^ $(GDAL_LIBS)
-
-*.o: global.h raster.h
-delineate_*.o: delineate_funcs.h
diff --git a/README.md b/README.md
index 47948bd..5fa8b07 100644
--- a/README.md
+++ b/README.md
@@ -13,12 +13,17 @@ Flow direction encoding in GeoTIFF:
* C compiler with [OpenMP](https://www.openmp.org/) support
* [GDAL](https://gdal.org/)
+* [Cmake](https://cmake.org/)
For Windows, use [MSYS2](https://www.msys2.org/) and [OSGeo4W](https://trac.osgeo.org/osgeo4w/) to install [GCC](https://gcc.gnu.org/) and [GDAL](https://gdal.org/), respectively.
## How to compile MESHED
+Use the standard cmake sequence:
```bash
+mkdir build
+cd build
+cmake -DCMAKE_BUILD_TYPE=RelWithDebInfo ..
make
```
diff --git a/delineate_funcs.h b/delineate_funcs.h
index 381483e..f3f621f 100644
--- a/delineate_funcs.h
+++ b/delineate_funcs.h
@@ -163,7 +163,7 @@ static void trace_up(struct raster_map *dir_map, int row, int col, int id,
/* use gcc -O2 or -O3 flags for tail-call optimization
* (-foptimize-sibling-calls) */
- trace_up(dir_map, row_next, col_next, id, up_stack);
+ MUST_TAIL return trace_up(dir_map, row_next, col_next, id, up_stack);
}
static void init_up_stack(struct cell_stack *up_stack)
diff --git a/global.h b/global.h
index 143f2ef..bdd4c76 100644
--- a/global.h
+++ b/global.h
@@ -83,4 +83,10 @@ struct hierarchy *analyze_hierarchy(struct raster_map *,
struct outlet_list *);
int write_hierarchy(const char *, struct hierarchy *);
+#if defined(__clang_major__) && __clang_major__ >= 13
+#define MUST_TAIL __attribute__((musttail))
+#else
+#define MUST_TAIL
+#endif
+
#endif
diff --git a/raster.c b/raster.c
index 035a791..cbcbf78 100644
--- a/raster.c
+++ b/raster.c
@@ -228,7 +228,7 @@ struct raster_map *read_raster(const char *path, int type, int get_stats)
#pragma omp parallel
{
#pragma omp single
- datasets = malloc(sizeof *datasets * omp_get_num_threads());
+ datasets = malloc(sizeof(GDALDatasetH*) * omp_get_num_threads());
datasets[omp_get_thread_num()] = GDALOpen(path, GA_ReadOnly);
}
@@ -246,7 +246,7 @@ struct raster_map *read_raster(const char *path, int type, int get_stats)
#pragma omp parallel
{
#pragma omp single
- bands = malloc(sizeof *bands * omp_get_num_threads());
+ bands = malloc(sizeof(GDALRasterBandH*) * omp_get_num_threads());
bands[omp_get_thread_num()] =
GDALGetRasterBand(datasets[omp_get_thread_num()], 1);
}
@@ -321,6 +321,9 @@ struct raster_map *read_raster(const char *path, int type, int get_stats)
if (error)
return NULL;
+ free(bands);
+ free(datasets);
+
return rast_map;
}
@@ -331,8 +334,6 @@ int write_raster(const char *path, struct raster_map *rast_map, int type)
GDALDatasetH dataset;
GDALRasterBandH band;
GDALDataType data_type, gdt_type;
- size_t row_size;
- int row;
if (!driver)
return 1;
@@ -342,9 +343,7 @@ int write_raster(const char *path, struct raster_map *rast_map, int type)
return 2;
if (rast_map->compress)
- options = CSLSetNameValue(options, "COMPRESS", "LZW");
-
- row_size = rast_map->ncols;
+ options = CSLSetNameValue(options, "COMPRESS", "ZSTD");
/* actual data size */
switch (rast_map->type) {
@@ -364,7 +363,6 @@ int write_raster(const char *path, struct raster_map *rast_map, int type)
data_type = GDT_Byte;
break;
}
- row_size *= GDALGetDataTypeSizeBytes(data_type);
/* requested data type */
gdt_type = data_type;
@@ -399,13 +397,13 @@ int write_raster(const char *path, struct raster_map *rast_map, int type)
band = GDALGetRasterBand(dataset, 1);
GDALSetRasterNoDataValue(band, rast_map->null_value);
- for (row = 0; row < rast_map->nrows; row++) {
- if (GDALRasterIO
- (band, GF_Write, 0, row, rast_map->ncols, 1,
- (char *)rast_map->cells.v + row * row_size, rast_map->ncols,
- 1, data_type, 0, 0) != CE_None)
- return 4;
- }
+ const int ret = GDALRasterIO(
+ band, GF_Write, 0, 0, rast_map->ncols, rast_map->nrows,
+ (char *)rast_map->cells.v, rast_map->ncols,
+ rast_map->nrows, data_type, 0, 0);
+
+ if (ret != CE_None)
+ return 4;
GDALClose(dataset);
diff --git a/timeval_diff.c b/timeval_diff.c
index cab3a0d..b32899e 100644
--- a/timeval_diff.c
+++ b/timeval_diff.c
@@ -1,4 +1,6 @@
#include
+#include "global.h"
+
#ifdef _MSC_VER
#include
#else