src/external/tng_io/src/tests/tng_io_read_pos_util.c

   1 /* This code is part of the tng binary trajectory format.
   2  *
   3  *  The high-level API of the TNG API is used where appropriate.
   4  *
   5  * Written by Magnus Lundborg
   6  * Copyright (c) 2012-2013, The GROMACS development team.
   7  * Check out http://www.gromacs.org for more information.
   8  *
   9  *
  10  * This program is free software; you can redistribute it and/or
  11  * modify it under the terms of the Revised BSD License.
  12  */
  13
  14 #ifdef USE_STD_INTTYPES_H
  15 #include <inttypes.h>
  16 #endif
  17
  18 #include <stdlib.h>
  19 #include <stdio.h>
  20 #include "../../include/tng_io.h"
  21
  22 int main(int argc, char **argv)
  23 {
  24     tng_trajectory_t traj;
  25     /* Assume that the data is stored as floats. The data is placed in 1-D
  26      * arrays */
  27     float *positions = 0, *box_shape = 0;
  28     int64_t n_particles, n_frames, tot_n_frames, stride_length, i, j;
  29     /* Set a default frame range */
  30     int first_frame = 0, last_frame = 5000, n_strides;
  31
  32     if(argc <= 1)
  33     {
  34         printf("No file specified\n");
  35         printf("Usage:\n");
  36         printf("tng_io_read_pos_util <tng_file> "
  37                "[first_frame = %d] [last_frame = %d]\n",
  38                first_frame, last_frame);
  39         exit(1);
  40     }
  41
  42     /* A reference must be passed to allocate memory */
  43     tng_util_trajectory_open(argv[1], 'r', &traj);
  44
  45     if(argc >= 3)
  46     {
  47         first_frame = strtol(argv[2], 0, 10);
  48         if(argc >= 4)
  49         {
  50             last_frame = strtol(argv[3], 0, 10);
  51         }
  52     }
  53
  54     if(tng_num_frames_get(traj, &tot_n_frames) != TNG_SUCCESS)
  55     {
  56         printf("Cannot determine the number of frames in the file\n");
  57         tng_util_trajectory_close(&traj);
  58         exit(1);
  59     }
  60
  61     if(tng_num_particles_get(traj, &n_particles) != TNG_SUCCESS)
  62     {
  63         printf("Cannot determine the number of particles in the file\n");
  64         tng_util_trajectory_close(&traj);
  65         exit(1);
  66     }
  67
  68     printf("%"PRId64" frames in file\n", tot_n_frames);
  69
  70     if(last_frame > tot_n_frames - 1)
  71     {
  72         last_frame = tot_n_frames - 1;
  73     }
  74
  75
  76     n_frames = last_frame - first_frame + 1;
  77
  78     if(tng_util_box_shape_read(traj, &box_shape, &stride_length) ==
  79        TNG_SUCCESS)
  80     {
  81         printf("Simulation box shape: ");
  82         for(i=0; i < 9; i++)
  83         {
  84             printf("%f ", box_shape[i]);
  85         }
  86         printf("\n");
  87     }
  88     else
  89     {
  90         printf("Simulation box shape not set in the file (or could not be read)\n");
  91     }
  92
  93     n_strides = (n_frames % stride_length) ? n_frames / stride_length + 1:
  94                 n_frames / stride_length;
  95
  96     /* Get the positions of all particles in the requested frame range.
  97      * The positions are stored in the positions array.
  98      * N.B. No proper error checks. */
  99     if(tng_util_pos_read_range(traj, 0, last_frame, &positions, &stride_length)
 100        == TNG_SUCCESS)
 101     {
 102         /* Print the positions of the wanted particle (zero based) */
 103         for(i=0; i < n_strides; i++)
 104         {
 105             printf("\nFrame %"PRId64":\n", first_frame + i*stride_length);
 106             for(j=0; j < n_particles; j++)
 107             {
 108                 printf("Atom nr: %"PRId64"", j);
 109                 printf("\t%f", positions[i*n_particles*3+j*3]);
 110                 printf("\t%f", positions[i*n_particles*3+j*3+1]);
 111                 printf("\t%f", positions[i*n_particles*3+j*3+2]);
 112                 printf("\n");
 113             }
 114         }
 115     }
 116     else
 117     {
 118         printf("Cannot read positions\n");
 119     }
 120
 121     /* Free memory */
 122     if(positions)
 123     {
 124         free(positions);
 125     }
 126     if(box_shape)
 127     {
 128         free(box_shape);
 129     }
 130     tng_util_trajectory_close(&traj);
 131
 132     return(0);
 133 }