c: iterate hdf5 datasets

CMakeLists.txt

cmake_minimum_required(VERSION 3.0)
project(hdf5_check_ds)

# output paths
set(CMAKE_ARCHIVE_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
set(CMAKE_LIBRARY_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/bin)

find_package(HDF5 REQUIRED COMPONENTS C HL)
include_directories(${HDF5_INCLUDE_DIRS})
link_libraries(${HDF5_C_LIBRARIES} ${HDF5_C_HL_LIBRARIES})
add_definitions(${HDF5_DEFINITIONS})
# THIS IS A TEMPORARY HACK TO OVERCOME A BUG IN SOME VERSIONS OF CMAKE
if(NOT (${HDF5_C_HL_LIBRARIES} MATCHES "libhdf5_hl"))
    list(GET HDF5_C_LIBRARIES 0 HDF5_LIB)
    get_filename_component(HDF5_LIBDIR ${HDF5_LIB} DIRECTORY)
    link_libraries("${HDF5_LIBDIR}/libhdf5_hl.so")
endif()

set(CMAKE_C_FLAGS "${cmake_c_flags} -std=gnu99 -g -O3 -march=native")

add_executable(hdf5_check_ds hdf5_check_ds.c)

hdf5_check_ds.c

#include <hdf5.h>
#include <hdf5_hl.h>
#include <stdlib.h>
#include <string.h>
#include <limits.h>

typedef struct unique_t {
  int *count;
  int max_id;
  int func_switch;
} unique_t;

int compare_ints(const void *a, const void *b) { return *(int *)a - *(int *)b; }

herr_t iterate_func(hid_t grp, const char *name, const H5L_info_t *info,
                    void *unique) {
  if (strstr(name, "Snap") != NULL) {
    // this is a snapshot group
    char ds_name[100];
    sprintf(ds_name, "%s/ForestID", name);

    hsize_t n_forests;
    H5LTget_dataset_info(grp, ds_name, &n_forests, NULL, NULL);

    int *forest_ids = NULL;
    if (n_forests > 0) {
      forest_ids = malloc(sizeof(int) * n_forests);
      H5LTread_dataset_int(grp, ds_name, forest_ids);

      if (((unique_t *)unique)->func_switch == 1) {
        int *count = ((unique_t *)unique)->count;
        for (int ii = 0; ii < (int)n_forests; ii++) count[forest_ids[ii] + 1]++;
      }

    } else if (((unique_t *)unique)->func_switch == 0) {
      n_forests = 2;
      forest_ids = malloc(sizeof(int) * n_forests);
      for (int ii = 0; ii < (int)n_forests; ii++) forest_ids[ii] = 0;
    }

    if (((unique_t *)unique)->func_switch == 0) {
      qsort(forest_ids, n_forests, sizeof(int), compare_ints);

      int *max_id = &(((unique_t *)unique)->max_id);
      if (forest_ids[n_forests - 1] > (int)(*max_id))
        *max_id = forest_ids[n_forests - 1];
    }

    if (forest_ids != NULL) free(forest_ids);
  }
  return 0;
}

int main() {
  const char fname[] = {
      "/lustre/projects/p124_astro/smutch/velociraptor/data/"
      "VELOCIraptor.tree.t4.unifiedhalotree.withforest.snap.hdf.data"};

  // open the file
  hid_t fd = H5Fopen(fname, H5F_ACC_RDONLY, H5P_DEFAULT);
  hid_t grp = H5Gopen(fd, "/", H5P_DEFAULT);

  // iterate through all members of the root group
  hsize_t idx = 0;
  unique_t unique;
  unique.func_switch = 0;
  unique.max_id = INT_MIN;
  H5Literate(grp, H5_INDEX_NAME, H5_ITER_INC, &idx, iterate_func, &unique);
  printf("max_id = %d\n", unique.max_id);

  int n_ids = unique.max_id + 2;
  unique.count = malloc(sizeof(int) * n_ids);
  for (int ii = 0; ii < n_ids; ii++) unique.count[ii] = 0;
  idx = 0;
  unique.func_switch = 1;
  H5Literate(grp, H5_INDEX_NAME, H5_ITER_INC, &idx, iterate_func, &unique);

  int n_unique = 0;
  for (int ii = 0; ii < n_ids; ii++)
    if (unique.count[ii] > 0) n_unique++;

  printf("n_unique = %d\n", n_unique);

  free(unique.count);

  // close the file
  H5Gclose(grp);
  H5Fclose(fd);

  return EXIT_SUCCESS;
}