Fix utility function for box shape reading.

Reading box shape using the high-level API should work.
Fixed bugs in data reading.
Updated utility usage examples.

2 files changed, 33 insertions, 17 deletions

diff --git a/src/tests/md_openmp_util.c b/src/tests/md_openmp_util.c
index 0196f98..339e1bb 100644
--- a/src/tests/md_openmp_util.c
+++ b/src/tests/md_openmp_util.c
@@ -118,6 +118,7 @@ int main ( int argc, char *argv[] )
 
     printf("\n");
     printf("  Initializing trajectory storage.\n");
+    /* Initialize the TNG trajectory */
 #ifdef EXAMPLE_FILES_DIR
     tng_util_trajectory_open(EXAMPLE_FILES_DIR "tng_md_out.tng", 'w', &traj);
 #else
@@ -127,6 +128,8 @@ int main ( int argc, char *argv[] )
 
 
     /* Set molecules data */
+    /* N.B. This is still not done using utility functions. The low-level API
+     * is used. */
     printf("  Creating molecules in trajectory.\n");
     tng_molecule_add(traj, "water", &molecule);
     tng_molecule_chain_add(traj, molecule, "W", &chain);
@@ -150,23 +153,11 @@ int main ( int argc, char *argv[] )
     for ( i = 0; i < nd; i++ )
     {
         box[i] = 10.0;
+        /* box_shape stores 9 values according to the TNG specs */
         box_shape[i*nd + i] = box[i];
     }
 
 
-    /* Add the box shape data block and write the file headers */
-//     if(tng_data_block_add(traj, TNG_TRAJ_BOX_SHAPE, "BOX SHAPE", TNG_DOUBLE_DATA,
-//                        TNG_NON_TRAJECTORY_BLOCK, 1, 9, 1, TNG_UNCOMPRESSED,
-//                        box_shape) == TNG_CRITICAL ||
-//                        tng_file_headers_write(traj, TNG_USE_HASH) == TNG_CRITICAL)
-//     {
-//         free(box_shape);
-//         tng_util_trajectory_close(&traj);
-//         printf("  Cannot write trajectory headers and box shape.\n");
-//         exit(1);
-//     }
-//     free(box_shape);
-
     printf ( "\n" );
     printf ( "  Initializing positions, velocities, and accelerations.\n" );
 /*
@@ -214,6 +205,15 @@ int main ( int argc, char *argv[] )
     step_print_index++;
     step_print = ( step_print_index * step_num ) / step_print_num;
 
+    /* The box shape does not change during the trajectory. */
+    if(tng_util_box_shape_write(traj, -1, box_shape) != TNG_SUCCESS)
+    {
+        printf("Error writing box shape. %s: %d\n",
+               __FILE__, __LINE__);
+        exit(1);
+    }
+
+    /* Set the output frequency of positions, velocities and forces */
     if(tng_util_pos_write_frequency_set(traj, step_save) != TNG_SUCCESS)
     {
         printf("Error setting writing frequency data. %s: %d\n",
@@ -233,6 +233,7 @@ int main ( int argc, char *argv[] )
         exit(1);
     }
 
+    /* Write the first frame of positions, velocities and forces */
     if(tng_util_pos_write(traj, 0, pos) != TNG_SUCCESS)
     {
         printf("Error adding data. %s: %d\n", __FILE__, __LINE__);
@@ -264,6 +265,7 @@ int main ( int argc, char *argv[] )
         }
         if(step % step_save == 0)
         {
+            /* Write positions, velocities and forces */
             if(tng_util_pos_write(traj, step, pos) != TNG_SUCCESS)
             {
                 printf("Error adding data. %s: %d\n", __FILE__, __LINE__);
@@ -294,9 +296,8 @@ int main ( int argc, char *argv[] )
     free ( force );
     free ( pos );
     free ( vel );
-/*
-    Terminate.
-*/
+
+    /* Close the TNG output. */
     tng_util_trajectory_close(&traj);
 
     printf ( "\n" );
diff --git a/src/tests/tng_io_read_pos_util.c b/src/tests/tng_io_read_pos_util.c
index 322f4a4..e7b00d3 100644
--- a/src/tests/tng_io_read_pos_util.c
+++ b/src/tests/tng_io_read_pos_util.c
@@ -24,7 +24,7 @@ int main(int argc, char **argv)
     tng_trajectory_t traj;
     /* Assume that the data is stored as floats. The data is placed in 1-D
      * arrays */
-    float *positions = 0, *velocities = 0, *forces = 0;
+    float *positions = 0, *box_shape = 0;
     int64_t n_particles, n_frames, tot_n_frames, stride_length, i, j;
     // Set a default frame range
     int64_t first_frame = 0, last_frame = 5000;
@@ -73,6 +73,21 @@ int main(int argc, char **argv)
         last_frame = tot_n_frames - 1;
     }
 
+    if(tng_util_box_shape_read(traj, &box_shape, &stride_length) ==
+        TNG_SUCCESS)
+    {
+        printf("Simulation box shape: ");
+        for(i=0; i < 9; i++)
+        {
+            printf("%f ", box_shape[i]);
+        }
+        printf("\n");
+    }
+    else
+    {
+        printf("Simulation box shape not set in the file (or could not be read)\n");
+    }
+
     n_frames = last_frame - first_frame + 1;