diff options
Diffstat (limited to 'src')
-rw-r--r-- | src/CMakeLists.txt | 2 | ||||
-rw-r--r-- | src/lib/CMakeLists.txt | 2 | ||||
-rw-r--r-- | src/lib/md5.c | 424 | ||||
-rw-r--r-- | src/lib/md5.h | 90 | ||||
-rw-r--r-- | src/lib/trg_io.c | 6934 | ||||
-rw-r--r-- | src/lib/trg_io.h | 760 | ||||
-rw-r--r-- | src/tests/CMakeLists.txt | 5 | ||||
-rw-r--r-- | src/tests/trg_io_testing.c | 341 |
8 files changed, 8558 insertions, 0 deletions
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt new file mode 100644 index 0000000..165c089 --- /dev/null +++ b/src/CMakeLists.txt @@ -0,0 +1,2 @@ +add_subdirectory(lib) +add_subdirectory(tests) diff --git a/src/lib/CMakeLists.txt b/src/lib/CMakeLists.txt new file mode 100644 index 0000000..c4c40bf --- /dev/null +++ b/src/lib/CMakeLists.txt @@ -0,0 +1,2 @@ +add_library(trg_io SHARED trg_io.c md5.c) + diff --git a/src/lib/md5.c b/src/lib/md5.c new file mode 100644 index 0000000..8bbde52 --- /dev/null +++ b/src/lib/md5.c @@ -0,0 +1,424 @@ +/* + * + * This source code is distributed with + * + * G R O M A C S + * + * GROningen MAchine for Chemical Simulations + * + * VERSION 4.5 + * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others. + * Copyright (c) 1991-2010, The GROMACS development team, + * check out http://www.gromacs.org for more information. + * + * Gromacs 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 2 + * of the License, or (at your option) any later version. + * + * If you want to redistribute modifications, please consider that + * scientific software is very special. Version control is crucial - + * bugs must be traceable. We will be happy to consider code for + * inclusion in the official distribution, but derived work must not + * be called official GROMACS. Details are found in the README & COPYING + * files - if they are missing, get the official version at www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the papers on the package - you can find them in the top README file. + * + * For more info, check our website at http://www.gromacs.org + * + * And Hey: + * GROningen Mixture of Alchemy and Childrens' Stories + */ +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#ifdef GMX_INTEGER_BIG_ENDIAN +#define ARCH_IS_BIG_ENDIAN 1 +#else +#define ARCH_IS_BIG_ENDIAN 0 +#endif + +/* + Copyright (C) 1999, 2000, 2002 Aladdin Enterprises. All rights reserved. + + This software is provided 'as-is', without any express or implied + warranty. In no event will the authors be held liable for any damages + arising from the use of this software. + + Permission is granted to anyone to use this software for any purpose, + including commercial applications, and to alter it and redistribute it + freely, subject to the following restrictions: + + 1. The origin of this software must not be misrepresented; you must not + claim that you wrote the original software. If you use this software + in a product, an acknowledgment in the product documentation would be + appreciated but is not required. + 2. Altered source versions must be plainly marked as such, and must not be + misrepresented as being the original software. + 3. This notice may not be removed or altered from any source distribution. + + L. Peter Deutsch + ghost@aladdin.com + + */ +/* + Independent implementation of MD5 (RFC 1321). + + This code implements the MD5 Algorithm defined in RFC 1321, whose + text is available at + http://www.ietf.org/rfc/rfc1321.txt + The code is derived from the text of the RFC, including the test suite + (section A.5) but excluding the rest of Appendix A. It does not include + any code or documentation that is identified in the RFC as being + copyrighted. + + The original and principal author of md5.c is L. Peter Deutsch + <ghost@aladdin.com>. Other authors are noted in the change history + that follows (in reverse chronological order): + + 2002-04-13 lpd Clarified derivation from RFC 1321; now handles byte order + either statically or dynamically; added missing #include <string.h> + in library. + 2002-03-11 lpd Corrected argument list for main(), and added int return + type, in test program and T value program. + 2002-02-21 lpd Added missing #include <stdio.h> in test program. + 2000-07-03 lpd Patched to eliminate warnings about "constant is + unsigned in ANSI C, signed in traditional"; made test program + self-checking. + 1999-11-04 lpd Edited comments slightly for automatic TOC extraction. + 1999-10-18 lpd Fixed typo in header comment (ansi2knr rather than md5). + 1999-05-03 lpd Original version. + */ + +#include "md5.h" +#include <string.h> + +#undef BYTE_ORDER /* 1 = big-endian, -1 = little-endian, 0 = unknown */ +#ifdef ARCH_IS_BIG_ENDIAN +# define BYTE_ORDER (ARCH_IS_BIG_ENDIAN ? 1 : -1) +#else +# define BYTE_ORDER 0 +#endif + +#define T_MASK ((md5_word_t)~0) +#define T1 /* 0xd76aa478 */ (T_MASK ^ 0x28955b87) +#define T2 /* 0xe8c7b756 */ (T_MASK ^ 0x173848a9) +#define T3 0x242070db +#define T4 /* 0xc1bdceee */ (T_MASK ^ 0x3e423111) +#define T5 /* 0xf57c0faf */ (T_MASK ^ 0x0a83f050) +#define T6 0x4787c62a +#define T7 /* 0xa8304613 */ (T_MASK ^ 0x57cfb9ec) +#define T8 /* 0xfd469501 */ (T_MASK ^ 0x02b96afe) +#define T9 0x698098d8 +#define T10 /* 0x8b44f7af */ (T_MASK ^ 0x74bb0850) +#define T11 /* 0xffff5bb1 */ (T_MASK ^ 0x0000a44e) +#define T12 /* 0x895cd7be */ (T_MASK ^ 0x76a32841) +#define T13 0x6b901122 +#define T14 /* 0xfd987193 */ (T_MASK ^ 0x02678e6c) +#define T15 /* 0xa679438e */ (T_MASK ^ 0x5986bc71) +#define T16 0x49b40821 +#define T17 /* 0xf61e2562 */ (T_MASK ^ 0x09e1da9d) +#define T18 /* 0xc040b340 */ (T_MASK ^ 0x3fbf4cbf) +#define T19 0x265e5a51 +#define T20 /* 0xe9b6c7aa */ (T_MASK ^ 0x16493855) +#define T21 /* 0xd62f105d */ (T_MASK ^ 0x29d0efa2) +#define T22 0x02441453 +#define T23 /* 0xd8a1e681 */ (T_MASK ^ 0x275e197e) +#define T24 /* 0xe7d3fbc8 */ (T_MASK ^ 0x182c0437) +#define T25 0x21e1cde6 +#define T26 /* 0xc33707d6 */ (T_MASK ^ 0x3cc8f829) +#define T27 /* 0xf4d50d87 */ (T_MASK ^ 0x0b2af278) +#define T28 0x455a14ed +#define T29 /* 0xa9e3e905 */ (T_MASK ^ 0x561c16fa) +#define T30 /* 0xfcefa3f8 */ (T_MASK ^ 0x03105c07) +#define T31 0x676f02d9 +#define T32 /* 0x8d2a4c8a */ (T_MASK ^ 0x72d5b375) +#define T33 /* 0xfffa3942 */ (T_MASK ^ 0x0005c6bd) +#define T34 /* 0x8771f681 */ (T_MASK ^ 0x788e097e) +#define T35 0x6d9d6122 +#define T36 /* 0xfde5380c */ (T_MASK ^ 0x021ac7f3) +#define T37 /* 0xa4beea44 */ (T_MASK ^ 0x5b4115bb) +#define T38 0x4bdecfa9 +#define T39 /* 0xf6bb4b60 */ (T_MASK ^ 0x0944b49f) +#define T40 /* 0xbebfbc70 */ (T_MASK ^ 0x4140438f) +#define T41 0x289b7ec6 +#define T42 /* 0xeaa127fa */ (T_MASK ^ 0x155ed805) +#define T43 /* 0xd4ef3085 */ (T_MASK ^ 0x2b10cf7a) +#define T44 0x04881d05 +#define T45 /* 0xd9d4d039 */ (T_MASK ^ 0x262b2fc6) +#define T46 /* 0xe6db99e5 */ (T_MASK ^ 0x1924661a) +#define T47 0x1fa27cf8 +#define T48 /* 0xc4ac5665 */ (T_MASK ^ 0x3b53a99a) +#define T49 /* 0xf4292244 */ (T_MASK ^ 0x0bd6ddbb) +#define T50 0x432aff97 +#define T51 /* 0xab9423a7 */ (T_MASK ^ 0x546bdc58) +#define T52 /* 0xfc93a039 */ (T_MASK ^ 0x036c5fc6) +#define T53 0x655b59c3 +#define T54 /* 0x8f0ccc92 */ (T_MASK ^ 0x70f3336d) +#define T55 /* 0xffeff47d */ (T_MASK ^ 0x00100b82) +#define T56 /* 0x85845dd1 */ (T_MASK ^ 0x7a7ba22e) +#define T57 0x6fa87e4f +#define T58 /* 0xfe2ce6e0 */ (T_MASK ^ 0x01d3191f) +#define T59 /* 0xa3014314 */ (T_MASK ^ 0x5cfebceb) +#define T60 0x4e0811a1 +#define T61 /* 0xf7537e82 */ (T_MASK ^ 0x08ac817d) +#define T62 /* 0xbd3af235 */ (T_MASK ^ 0x42c50dca) +#define T63 0x2ad7d2bb +#define T64 /* 0xeb86d391 */ (T_MASK ^ 0x14792c6e) + + +static inline void +md5_process(md5_state_t *pms, const md5_byte_t *data /*[64]*/) +{ + md5_word_t + a = pms->abcd[0], b = pms->abcd[1], + c = pms->abcd[2], d = pms->abcd[3]; + md5_word_t t; +#if BYTE_ORDER > 0 + /* Define storage only for big-endian CPUs. */ + md5_word_t X[16]; +#else + /* Define storage for little-endian or both types of CPUs. */ + md5_word_t xbuf[16]; + /* cppcheck-suppress unassignedVariable */ + const md5_word_t *X; +#endif + + { +#if BYTE_ORDER == 0 + /* + * Determine dynamically whether this is a big-endian or + * little-endian machine, since we can use a more efficient + * algorithm on the latter. + */ + static const int w = 1; + + if (*((const md5_byte_t *)&w)) /* dynamic little-endian */ +#endif +#if BYTE_ORDER <= 0 /* little-endian */ + { + /* + * On little-endian machines, we can process properly aligned + * data without copying it. + */ + if (!((data - (const md5_byte_t *)0) & 3)) { + /* data are properly aligned */ + X = (const md5_word_t *)data; + } else { + /* not aligned */ + memcpy(xbuf, data, 64); + X = xbuf; + } + } +#endif +#if BYTE_ORDER == 0 + else /* dynamic big-endian */ +#endif +#if BYTE_ORDER >= 0 /* big-endian */ + { + /* + * On big-endian machines, we must arrange the bytes in the + * right order. + */ + const md5_byte_t *xp = data; + int i; + +# if BYTE_ORDER == 0 + X = xbuf; /* (dynamic only) */ +# else +# define xbuf X /* (static only) */ +# endif + for (i = 0; i < 16; ++i, xp += 4) + xbuf[i] = xp[0] + (xp[1] << 8) + (xp[2] << 16) + (xp[3] << 24); + } +#endif + } + +#define ROTATE_LEFT(x, n) (((x) << (n)) | ((x) >> (32 - (n)))) + + /* Round 1. */ + /* Let [abcd k s i] denote the operation + a = b + ((a + F(b,c,d) + X[k] + T[i]) <<< s). */ +#define F(x, y, z) (((x) & (y)) | (~(x) & (z))) +#define SET(a, b, c, d, k, s, Ti)\ + t = a + F(b,c,d) + X[k] + Ti;\ + a = ROTATE_LEFT(t, s) + b + /* Do the following 16 operations. */ + SET(a, b, c, d, 0, 7, T1); + SET(d, a, b, c, 1, 12, T2); + SET(c, d, a, b, 2, 17, T3); + SET(b, c, d, a, 3, 22, T4); + SET(a, b, c, d, 4, 7, T5); + SET(d, a, b, c, 5, 12, T6); + SET(c, d, a, b, 6, 17, T7); + SET(b, c, d, a, 7, 22, T8); + SET(a, b, c, d, 8, 7, T9); + SET(d, a, b, c, 9, 12, T10); + SET(c, d, a, b, 10, 17, T11); + SET(b, c, d, a, 11, 22, T12); + SET(a, b, c, d, 12, 7, T13); + SET(d, a, b, c, 13, 12, T14); + SET(c, d, a, b, 14, 17, T15); + SET(b, c, d, a, 15, 22, T16); +#undef SET + + /* Round 2. */ + /* Let [abcd k s i] denote the operation + a = b + ((a + G(b,c,d) + X[k] + T[i]) <<< s). */ +#define G(x, y, z) (((x) & (z)) | ((y) & ~(z))) +#define SET(a, b, c, d, k, s, Ti)\ + t = a + G(b,c,d) + X[k] + Ti;\ + a = ROTATE_LEFT(t, s) + b + /* Do the following 16 operations. */ + SET(a, b, c, d, 1, 5, T17); + SET(d, a, b, c, 6, 9, T18); + SET(c, d, a, b, 11, 14, T19); + SET(b, c, d, a, 0, 20, T20); + SET(a, b, c, d, 5, 5, T21); + SET(d, a, b, c, 10, 9, T22); + SET(c, d, a, b, 15, 14, T23); + SET(b, c, d, a, 4, 20, T24); + SET(a, b, c, d, 9, 5, T25); + SET(d, a, b, c, 14, 9, T26); + SET(c, d, a, b, 3, 14, T27); + SET(b, c, d, a, 8, 20, T28); + SET(a, b, c, d, 13, 5, T29); + SET(d, a, b, c, 2, 9, T30); + SET(c, d, a, b, 7, 14, T31); + SET(b, c, d, a, 12, 20, T32); +#undef SET + + /* Round 3. */ + /* Let [abcd k s t] denote the operation + a = b + ((a + H(b,c,d) + X[k] + T[i]) <<< s). */ +#define H(x, y, z) ((x) ^ (y) ^ (z)) +#define SET(a, b, c, d, k, s, Ti)\ + t = a + H(b,c,d) + X[k] + Ti;\ + a = ROTATE_LEFT(t, s) + b + /* Do the following 16 operations. */ + SET(a, b, c, d, 5, 4, T33); + SET(d, a, b, c, 8, 11, T34); + SET(c, d, a, b, 11, 16, T35); + SET(b, c, d, a, 14, 23, T36); + SET(a, b, c, d, 1, 4, T37); + SET(d, a, b, c, 4, 11, T38); + SET(c, d, a, b, 7, 16, T39); + SET(b, c, d, a, 10, 23, T40); + SET(a, b, c, d, 13, 4, T41); + SET(d, a, b, c, 0, 11, T42); + SET(c, d, a, b, 3, 16, T43); + SET(b, c, d, a, 6, 23, T44); + SET(a, b, c, d, 9, 4, T45); + SET(d, a, b, c, 12, 11, T46); + SET(c, d, a, b, 15, 16, T47); + SET(b, c, d, a, 2, 23, T48); +#undef SET + + /* Round 4. */ + /* Let [abcd k s t] denote the operation + a = b + ((a + I(b,c,d) + X[k] + T[i]) <<< s). */ +#define I(x, y, z) ((y) ^ ((x) | ~(z))) +#define SET(a, b, c, d, k, s, Ti)\ + t = a + I(b,c,d) + X[k] + Ti;\ + a = ROTATE_LEFT(t, s) + b + /* Do the following 16 operations. */ + SET(a, b, c, d, 0, 6, T49); + SET(d, a, b, c, 7, 10, T50); + SET(c, d, a, b, 14, 15, T51); + SET(b, c, d, a, 5, 21, T52); + SET(a, b, c, d, 12, 6, T53); + SET(d, a, b, c, 3, 10, T54); + SET(c, d, a, b, 10, 15, T55); + SET(b, c, d, a, 1, 21, T56); + SET(a, b, c, d, 8, 6, T57); + SET(d, a, b, c, 15, 10, T58); + SET(c, d, a, b, 6, 15, T59); + SET(b, c, d, a, 13, 21, T60); + SET(a, b, c, d, 4, 6, T61); + SET(d, a, b, c, 11, 10, T62); + SET(c, d, a, b, 2, 15, T63); + SET(b, c, d, a, 9, 21, T64); +#undef SET + + /* Then perform the following additions. (That is increment each + of the four registers by the value it had before this block + was started.) */ + pms->abcd[0] += a; + pms->abcd[1] += b; + pms->abcd[2] += c; + pms->abcd[3] += d; +} + +void +md5_init(md5_state_t *pms) +{ + pms->count[0] = pms->count[1] = 0; + pms->abcd[0] = 0x67452301; + pms->abcd[1] = /*0xefcdab89*/ T_MASK ^ 0x10325476; + pms->abcd[2] = /*0x98badcfe*/ T_MASK ^ 0x67452301; + pms->abcd[3] = 0x10325476; +} + +void +md5_append(md5_state_t *pms, const md5_byte_t *data, int nbytes) +{ + const md5_byte_t *p = data; + int left = nbytes; + int offset = (pms->count[0] >> 3) & 63; + md5_word_t nbits = (md5_word_t)(nbytes << 3); + + if (nbytes <= 0) + return; + + /* Update the message length. */ + pms->count[1] += nbytes >> 29; + pms->count[0] += nbits; + if (pms->count[0] < nbits) + pms->count[1]++; + + /* Process an initial partial block. */ + if (offset) { + int copy = (offset + nbytes > 64 ? 64 - offset : nbytes); + + memcpy(pms->buf + offset, p, copy); + if (offset + copy < 64) + return; + p += copy; + left -= copy; + md5_process(pms, pms->buf); + } + + /* Process full blocks. */ + for (; left >= 64; p += 64, left -= 64) + md5_process(pms, p); + + /* Process a final partial block. */ + if (left) + memcpy(pms->buf, p, left); +} + +void +md5_finish(md5_state_t *pms, md5_byte_t digest[16]) +{ + static const md5_byte_t pad[64] = { + 0x80, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 + }; + md5_byte_t data[8]; + int i; + + /* Save the length before padding. */ + for (i = 0; i < 8; ++i) + data[i] = (md5_byte_t)(pms->count[i >> 2] >> ((i & 3) << 3)); + /* Pad to 56 bytes mod 64. */ + md5_append(pms, pad, ((55 - (pms->count[0] >> 3)) & 63) + 1); + /* Append the length. */ + md5_append(pms, data, 8); + for (i = 0; i < 16; ++i) + digest[i] = (md5_byte_t)(pms->abcd[i >> 2] >> ((i & 3) << 3)); +} diff --git a/src/lib/md5.h b/src/lib/md5.h new file mode 100644 index 0000000..15719fd --- /dev/null +++ b/src/lib/md5.h @@ -0,0 +1,90 @@ +/* + Copyright (C) 1999, 2002 Aladdin Enterprises. All rights reserved. + + This software is provided 'as-is', without any express or implied + warranty. In no event will the authors be held liable for any damages + arising from the use of this software. + + Permission is granted to anyone to use this software for any purpose, + including commercial applications, and to alter it and redistribute it + freely, subject to the following restrictions: + + 1. The origin of this software must not be misrepresented; you must not + claim that you wrote the original software. If you use this software + in a product, an acknowledgment in the product documentation would be + appreciated but is not required. + 2. Altered source versions must be plainly marked as such, and must not be + misrepresented as being the original software. + 3. This notice may not be removed or altered from any source distribution. + + L. Peter Deutsch + ghost@aladdin.com + + */ +/* + Independent implementation of MD5 (RFC 1321). + + This code implements the MD5 Algorithm defined in RFC 1321, whose + text is available at + http://www.ietf.org/rfc/rfc1321.txt + The code is derived from the text of the RFC, including the test suite + (section A.5) but excluding the rest of Appendix A. It does not include + any code or documentation that is identified in the RFC as being + copyrighted. + + The original and principal author of md5.h is L. Peter Deutsch + <ghost@aladdin.com>. Other authors are noted in the change history + that follows (in reverse chronological order): + + 2002-04-13 lpd Removed support for non-ANSI compilers; removed + references to Ghostscript; clarified derivation from RFC 1321; + now handles byte order either statically or dynamically. + 1999-11-04 lpd Edited comments slightly for automatic TOC extraction. + 1999-10-18 lpd Fixed typo in header comment (ansi2knr rather than md5); + added conditionalization for C++ compilation from Martin + Purschke <purschke@bnl.gov>. + 1999-05-03 lpd Original version. + */ + +#ifndef md5_INCLUDED +# define md5_INCLUDED + +/* + * This package supports both compile-time and run-time determination of CPU + * byte order. If ARCH_IS_BIG_ENDIAN is defined as 0, the code will be + * compiled to run only on little-endian CPUs; if ARCH_IS_BIG_ENDIAN is + * defined as non-zero, the code will be compiled to run only on big-endian + * CPUs; if ARCH_IS_BIG_ENDIAN is not defined, the code will be compiled to + * run on either big- or little-endian CPUs, but will run slightly less + * efficiently on either one than if ARCH_IS_BIG_ENDIAN is defined. + */ + +typedef unsigned char md5_byte_t; /* 8-bit byte */ +typedef unsigned int md5_word_t; /* 32-bit word */ + +/* Define the state of the MD5 Algorithm. */ +typedef struct md5_state_s { + md5_word_t count[2]; /* message length in bits, lsw first */ + md5_word_t abcd[4]; /* digest buffer */ + md5_byte_t buf[64]; /* accumulate block */ +} md5_state_t; + +#ifdef __cplusplus +extern "C" +{ +#endif + +/* Initialize the algorithm. */ +void md5_init(md5_state_t *pms); + +/* Append a string to the message. */ +void md5_append(md5_state_t *pms, const md5_byte_t *data, int nbytes); + +/* Finish the message and return the digest. */ +void md5_finish(md5_state_t *pms, md5_byte_t digest[16]); + +#ifdef __cplusplus +} /* end extern "C" */ +#endif + +#endif /* md5_INCLUDED */ diff --git a/src/lib/trg_io.c b/src/lib/trg_io.c new file mode 100644 index 0000000..c05d078 --- /dev/null +++ b/src/lib/trg_io.c @@ -0,0 +1,6934 @@ +/* This code is part of the trg binary trajectory format. + * + * VERSION 1.0 + * + * Written by Magnus Lundborg + * Copyright (c) 2012, The GROMACS development team. + * check out http://www.gromacs.org for more information. + * + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + */ + + +/* TODO: Make sure UTF-8 works properly */ + +#include <inttypes.h> +#include <limits.h> +#include <stdlib.h> +#include <string.h> +#include <time.h> +#include <unistd.h> + +#include "trg_io.h" +#include "md5.h" + +#ifndef LOGIN_NAME_MAX +#define LOGIN_NAME_MAX 33 +#endif +#ifndef HOST_NAME_MAX +#define HOST_NAME_MAX 257 +#endif + + +/* This function swaps the byte order of a 32 bit numerical variable. + It does not only work with integer, but e.g. floats need casting */ +static inline trg_function_status trg_swap_byte_order_32 + (const struct trg_trajectory *trg_data, uint32_t *v) +{ + switch(trg_data->endianness_32) + { + case TRG_LITTLE_ENDIAN_32: // Byte order is reversed. + *v = (*v << 24) | // Move first byte to end + ((*v & 0x00FF0000) >> 8) | // Move 2nd byte to pos 3 + ((*v & 0x0000FF00) << 8) | // Move 3rd byte to pos 2 + (*v >> 24); // Move last byte to first + + return(TRG_SUCCESS); + + case TRG_BYTE_PAIR_SWAP_32: // byte pair swap + *v = ((*v & 0xFFFF0000) >> 16) | + ((*v & 0x0000FFFF) << 16); + + return(TRG_SUCCESS); + + case TRG_BIG_ENDIAN_32: // Already correct + return(TRG_SUCCESS); + + default: + return(TRG_FAILURE); + } +} + +/* This function swaps the byte order of a 64 bit numerical variable. + It does not only work with integer, but e.g. floats need casting + The byte order swapping routine can convert five different byte + orders to big endian. */ +static inline trg_function_status trg_swap_byte_order_64 + (const struct trg_trajectory *trg_data, uint64_t *v) +{ + switch(trg_data->endianness_64) + { + case TRG_LITTLE_ENDIAN_64: // Byte order is reversed. + *v = (*v >> 56) | // Move first byte to end + ((*v & 0x00FF000000000000) >> 40) | // Move 2nd byte to pos 7 + ((*v & 0x0000FF0000000000) >> 24) | // Move 3rd byte to pos 6 + ((*v & 0x000000FF00000000) >> 8 ) | // Move 4th byte to pos 5 + ((*v & 0x00000000FF000000) << 8 ) | // Move 5th byte to pos 4 + ((*v & 0x0000000000FF0000) << 24) | // Move 6th byte to pos 3 + ((*v & 0x000000000000FF00) << 40) | // Move 7th byte to pos 2 + (*v << 56); // Move last byte to first + + return(TRG_SUCCESS); + + case TRG_QUAD_SWAP_64: // Byte quad swap + *v = ((*v & 0xFFFFFFFF00000000) >> 32) | + ((*v & 0x00000000FFFFFFFF) << 32); + + return(TRG_SUCCESS); + + case TRG_BYTE_PAIR_SWAP_64: // Byte pair swap + *v = ((*v & 0xFFFF0000FFFF0000) >> 16) | + ((*v & 0x0000FFFF0000FFFF) << 16); + + return(TRG_SUCCESS); + + case TRG_BYTE_SWAP_64: // Byte swap + *v = ((*v & 0xFF00FF00FF00FF00) >> 8) | + ((*v & 0x00FF00FF00FF00FF) << 8); + + return(TRG_SUCCESS); + + case TRG_BIG_ENDIAN_64: // Already correct + return(TRG_SUCCESS); + + default: + return(TRG_FAILURE); + } +} + +/* Generate the md5 hash of a block. + The hash is created based on the actual block contents. */ +static trg_function_status trg_generate_block_hash(struct trg_gen_block *block) +{ + md5_state_t md5_state; + + md5_init(&md5_state); + md5_append(&md5_state, (md5_byte_t *)block->block_contents, + block->block_contents_size); + md5_finish(&md5_state, (md5_byte_t *)block->hash); + + return(TRG_SUCCESS); +} + +/* Compare the current block hash (e.g. read from file) with the hash + calculated from the current contents. + If the current hash is all zeros skip the comparison. + If the hashes match results is set to TRUE, otherwise it is set to + FALSE. */ +static trg_function_status trg_verify_hash_match(struct trg_gen_block *block, + trg_bool *results) +{ + md5_state_t md5_state; + char hash[TRG_HASH_LEN]; + + if(strcmp(block->hash,"\0\0\0\0\0\0\0\0\0\0\0\0\0\0\0\0") == 0) + { + *results = TRUE; + return(TRG_SUCCESS); + } + md5_init(&md5_state); + md5_append(&md5_state, (md5_byte_t *)block->block_contents, + block->block_contents_size); + md5_finish(&md5_state, (md5_byte_t *)hash); + + if(strncmp(block->hash, hash, 16) != 0) + { + *results = FALSE; + } + else + { + *results = TRUE; + } + + return(TRG_SUCCESS); +} + + +static trg_function_status trg_init_input_file(struct trg_trajectory *trg_data, + const trg_bool update_read_pos) +{ + if(!trg_data->input_file) + { + if(!trg_data->input_file_path) + { + printf("No file specified for reading. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + trg_data->input_file = fopen(trg_data->input_file_path, "r"); + if(!trg_data->input_file) + { + printf("Cannot open file %s. %s: %d\n", + trg_data->input_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(fseek(trg_data->input_file, trg_data->input_file_pos, + SEEK_SET) != 0) + { + printf("Cannot specify position in file %s. %s: %d\n", + trg_data->input_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + else if(update_read_pos && fseek(trg_data->input_file, + trg_data->input_file_pos, SEEK_SET) != 0) + { + printf("Cannot specify position in file %s. %s: %d\n", + trg_data->input_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); +} + +static trg_function_status trg_init_output_file + (struct trg_trajectory *trg_data, + const trg_bool update_write_pos) +{ + if(!trg_data->output_file) + { + if(!trg_data->output_file_path) + { + printf("No file specified for writing. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + if(trg_data->output_file_pos <= 0) + { + trg_data->output_file = fopen(trg_data->output_file_path, "w+"); + } + else + { + trg_data->output_file = fopen(trg_data->output_file_path, "a+"); + } + + if(!trg_data->output_file) + { + printf("Cannot open file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(fseek(trg_data->output_file, 0, SEEK_SET) != 0) + { + printf("Cannot specify position in file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + else if(update_write_pos && fseek(trg_data->output_file, 0, + SEEK_SET) != 0) + { + printf("Cannot specify position in file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_block_header + (struct trg_trajectory *trg_data, struct trg_gen_block *block) +{ + int len, offset = 0; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + /* First read the header size to be able to read the whole header. */ + if(fread(&block->header_contents_size, sizeof(block->header_contents_size), + 1, trg_data->input_file) == 0) + { + printf("Cannot read header size. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &block->header_contents_size) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + /* Move the reading position to the beginning of the header. */ + fseek(trg_data->input_file, -sizeof(block->header_contents_size), + SEEK_CUR); + + /* If there is already memory allocated for the contents free it (we do not + * know if the size is correct). */ + if(block->header_contents) + { + free(block->header_contents); + } + + block->header_contents = (char *) malloc(block->header_contents_size); + if(!block->header_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->header_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole header into header_contents. This way it can be saved + * even if it cannot be interpreted + * for one reason or another. */ + if(fread(block->header_contents, block->header_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read header. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* The header contents size has already been read. Skip ahead. */ + offset = sizeof(block->header_contents_size); + + + /* Copy the respective parameters from the header contents block */ + memcpy(&block->block_contents_size, block->header_contents+offset, + sizeof(block->block_contents_size)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &block->block_contents_size) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->block_contents_size); + + memcpy(&block->id, block->header_contents+offset, sizeof(block->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &block->id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->id); + + memcpy(block->hash, block->header_contents+offset, TRG_HASH_LEN); + offset += TRG_HASH_LEN; + + if(block->name && strcmp(block->name, block->header_contents+offset) != 0) + { + free(block->name); + block->name = 0; + } + len = min(strlen(block->header_contents+offset) + 1, TRG_MAX_STR_LEN); + if(!block->name) + { + block->name = (char *) malloc(len); + if(!block->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strncpy(block->name, block->header_contents+offset, len); + } + offset += len; + + memcpy(&block->block_version, block->header_contents+offset, + sizeof(block->block_version)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &block->block_version) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->block_version); + + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_block_header + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, write_mode mode) +{ + int name_len, offset = 0; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + printf("Cannot initialise destination file. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + + /* Just dump the full header contents to file */ + if(mode == TRG_COPY_EXISTING) + { + if(!block->header_contents) + { + printf("No contents to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->header_contents, block->header_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all header data. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + if(!block->name) + { + block->name = (char *) malloc(1); + if(!block->name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + block->name[0] = 0; + } + + name_len = min(strlen(block->name) + 1, TRG_MAX_STR_LEN); + + trg_generate_block_hash(block); + + /* Calculate the size of the header to write */ + block->header_contents_size = sizeof(block->header_contents_size) + + sizeof(block->block_contents_size) + + sizeof(block->id) + + sizeof(block->block_version) + + TRG_HASH_LEN + + name_len; + + if(block->header_contents) + { + free(block->header_contents); + } + + block->header_contents = (char *) malloc(block->header_contents_size); + if(!block->header_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->header_contents_size, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + /* First copy all data into the header_contents block and finally write + * the whole block at once. */ + memcpy(block->header_contents, &block->header_contents_size, + sizeof(block->header_contents_size)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->header_contents)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->header_contents_size); + + memcpy(block->header_contents+offset, &block->block_contents_size, + sizeof(block->block_contents_size)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->header_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->block_contents_size); + + memcpy(block->header_contents+offset, &block->id, sizeof(block->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->header_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->id); + + memcpy(block->header_contents+offset, block->hash, TRG_HASH_LEN); + offset += TRG_HASH_LEN; + + strncpy(block->header_contents+offset, block->name, name_len); + offset += name_len; + + memcpy(block->header_contents+offset, &block->block_version, + sizeof(block->block_version)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->header_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block->block_version); + + if(fwrite(block->header_contents, block->header_contents_size, + 1, trg_data->output_file) != 1) + { + printf("Could not write all header data. %s: %d\n", __FILE__, __LINE__); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); +} + + +static trg_function_status trg_read_general_info_block + (struct trg_trajectory *trg_data, struct trg_gen_block *block) +{ + int len, offset = 0; + trg_bool same_hash; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to disk + * even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("General info block contents corrupt. Hashes do not match. " + "%s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + len = min(strlen(block->block_contents) + 1, TRG_MAX_STR_LEN); + trg_data->program_name = (char *) malloc(len); + if(!trg_data->program_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strncpy(trg_data->program_name, block->block_contents, len); + offset += len; + + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + trg_data->forcefield_name = (char *) malloc(len); + if(!trg_data->forcefield_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strncpy(trg_data->forcefield_name, block->block_contents+offset, len); + offset += len; + + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + trg_data->user_name = (char *) malloc(len); + if(!trg_data->user_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strncpy(trg_data->user_name, block->block_contents+offset, len); + offset += len; + + memcpy(&trg_data->time, block->block_contents+offset, + sizeof(trg_data->time)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &trg_data->time) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->time); + + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + trg_data->computer_name = (char *) malloc(len); + if(!trg_data->computer_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(trg_data->computer_name, block->block_contents+offset, len); + offset += len; + + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + trg_data->pgp_signature = (char *) malloc(len); + if(!trg_data->pgp_signature) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(trg_data->pgp_signature, block->block_contents+offset, len); + offset += len; + + memcpy(&trg_data->var_num_atoms_flag, block->block_contents+offset, + sizeof(trg_data->var_num_atoms_flag)); + offset += sizeof(trg_data->var_num_atoms_flag); + + memcpy(&trg_data->frame_set_n_frames, block->block_contents+offset, + sizeof(trg_data->frame_set_n_frames)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &trg_data->frame_set_n_frames) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->frame_set_n_frames); + + memcpy(&trg_data->first_trajectory_frame_set_input_file_pos, + block->block_contents+offset, + sizeof(trg_data->first_trajectory_frame_set_input_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &trg_data->first_trajectory_frame_set_input_file_pos) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->first_trajectory_frame_set_input_file_pos); + trg_data->current_trajectory_frame_set.next_frame_set_file_pos = + trg_data->first_trajectory_frame_set_input_file_pos; + + + memcpy(&trg_data->last_trajectory_frame_set_input_file_pos, + block->block_contents+offset, + sizeof(trg_data->last_trajectory_frame_set_input_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &trg_data->last_trajectory_frame_set_input_file_pos) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->last_trajectory_frame_set_input_file_pos); + + memcpy(&trg_data->stride_length, block->block_contents+offset, + sizeof(trg_data->stride_length)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &trg_data->stride_length) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_general_info_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + write_mode mode) +{ + int program_name_len, forcefield_name_len, user_name_len; + int computer_name_len, pgp_signature_len; + int offset = 0; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + if(!trg_data->program_name) + { + trg_data->program_name = (char *) malloc(1); + if(!trg_data->program_name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + trg_data->program_name[0] = 0; + } + if(!trg_data->forcefield_name) + { + trg_data->forcefield_name = (char *) malloc(1); + if(!trg_data->forcefield_name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + trg_data->forcefield_name[0] = 0; + } + if(!trg_data->user_name) + { + trg_data->user_name = (char *) malloc(1); + if(!trg_data->user_name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + trg_data->user_name[0] = 0; + } + if(!trg_data->computer_name) + { + trg_data->computer_name = (char *) malloc(1); + if(!trg_data->computer_name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + trg_data->computer_name[0] = 0; + } + if(!trg_data->pgp_signature) + { + trg_data->pgp_signature = (char *) malloc(1); + if(!trg_data->pgp_signature) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + trg_data->pgp_signature[0] = 0; + } + + program_name_len = min(strlen(trg_data->program_name) + 1, + TRG_MAX_STR_LEN); + forcefield_name_len = min(strlen(trg_data->forcefield_name) + 1, + TRG_MAX_STR_LEN); + user_name_len = min(strlen(trg_data->user_name) + 1, + TRG_MAX_STR_LEN); + computer_name_len = min(strlen(trg_data->computer_name) + 1, + TRG_MAX_STR_LEN); + pgp_signature_len = min(strlen(trg_data->pgp_signature) + 1, + TRG_MAX_STR_LEN); + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", + __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + block->block_contents_size = sizeof(trg_data->time) + + sizeof(trg_data->var_num_atoms_flag) + + sizeof(trg_data->frame_set_n_frames) + + sizeof(trg_data->first_trajectory_frame_set_input_file_pos) + + sizeof(trg_data->last_trajectory_frame_set_input_file_pos) + + sizeof(trg_data->stride_length) + + program_name_len + + forcefield_name_len + + user_name_len + + computer_name_len + + pgp_signature_len; + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + strncpy(block->block_contents, trg_data->program_name, program_name_len); + offset += program_name_len; + + strncpy(block->block_contents+offset, trg_data->forcefield_name, + forcefield_name_len); + offset += forcefield_name_len; + + strncpy(block->block_contents+offset, trg_data->user_name, user_name_len); + offset += user_name_len; + + memcpy(block->block_contents+offset, &trg_data->time, + sizeof(trg_data->time)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->time); + + strncpy(block->block_contents+offset, trg_data->computer_name, + computer_name_len); + offset += computer_name_len; + + strncpy(block->block_contents+offset, trg_data->pgp_signature, + pgp_signature_len); + offset += pgp_signature_len; + + memcpy(block->block_contents+offset, &trg_data->var_num_atoms_flag, + sizeof(trg_data->var_num_atoms_flag)); + offset += sizeof(trg_data->var_num_atoms_flag); + + memcpy(block->block_contents+offset, &trg_data->frame_set_n_frames, + sizeof(trg_data->frame_set_n_frames)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->frame_set_n_frames); + + memcpy(block->block_contents+offset, + &trg_data->first_trajectory_frame_set_input_file_pos, + sizeof(trg_data->first_trajectory_frame_set_input_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->first_trajectory_frame_set_input_file_pos); + + memcpy(block->block_contents+offset, + &trg_data->last_trajectory_frame_set_input_file_pos, + sizeof(trg_data->last_trajectory_frame_set_input_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->last_trajectory_frame_set_input_file_pos); + + + memcpy(block->block_contents+offset, &trg_data->stride_length, + sizeof(trg_data->stride_length)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +/* FIXME: Update this according to the new specs */ +static trg_function_status trg_read_molecules_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block) +{ + int i, j, k, l, len, offset = 0; + struct trg_molecule *molecule; + struct trg_chain *chain; + struct trg_residue *residue; + struct trg_atom *atom; + struct trg_bond *bond; + trg_bool same_hash; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to disk + * even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("Molecules block contents corrupt. Hashes do not match. " + "%s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + memcpy(&trg_data->n_molecules, block->block_contents, + sizeof(trg_data->n_molecules)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &trg_data->n_molecules) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->n_molecules); + + if(trg_data->molecules) + { + free(trg_data->molecules); + } + + trg_data->n_particles = 0; + + trg_data->molecules = (struct trg_molecule *) + malloc(trg_data->n_molecules * + sizeof(struct trg_molecule)); + if(!trg_data->molecules) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + trg_data->n_molecules * sizeof(struct trg_molecule), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!trg_data->var_num_atoms_flag) + { + if(trg_data->molecule_cnt_list) + { + free(trg_data->molecule_cnt_list); + } + trg_data->molecule_cnt_list = (uint64_t *) malloc(sizeof(uint64_t) * + trg_data->n_molecules); + if(!trg_data->molecule_cnt_list) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + trg_data->n_molecules * sizeof(struct trg_molecule), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + } + + + for(i=0; i < trg_data->n_molecules; i++) + { + molecule = &trg_data->molecules[i]; + + memcpy(&molecule->id, block->block_contents+offset, + sizeof(molecule->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->id); + +// printf("Read id: %"PRId64" offset: %d\n", molecule->id, offset); + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + molecule->name = (char *) malloc(len); + strncpy(molecule->name, block->block_contents+offset, len); + offset += len; + + memcpy(&molecule->quaternary_str, block->block_contents+offset, + sizeof(molecule->quaternary_str)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->quaternary_str) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->quaternary_str); + + if(!trg_data->var_num_atoms_flag) + { + memcpy(&trg_data->molecule_cnt_list[i], + block->block_contents+offset, + sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &trg_data->molecule_cnt_list[i]) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + } + + + memcpy(&molecule->n_chains, block->block_contents+offset, + sizeof(molecule->n_chains)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->n_chains) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_chains); + + memcpy(&molecule->n_residues, block->block_contents+offset, + sizeof(molecule->n_residues)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->n_residues) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_residues); + + memcpy(&molecule->n_atoms, block->block_contents+offset, + sizeof(molecule->n_atoms)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->n_atoms) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_atoms); + + trg_data->n_particles += molecule->n_atoms * + trg_data->molecule_cnt_list[i]; + + molecule->chains = (struct trg_chain *) malloc(molecule->n_chains * + sizeof(struct trg_chain)); + if(!molecule->chains) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + molecule->n_chains * sizeof(struct trg_chain), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + chain = molecule->chains; + + molecule->residues = (struct trg_residue *) + malloc(molecule->n_residues * + sizeof(struct trg_residue)); + if(!molecule->residues) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + molecule->n_residues * sizeof(struct trg_residue), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + residue = molecule->residues; + + molecule->atoms = (struct trg_atom *) malloc(molecule->n_atoms * + sizeof(struct trg_atom)); + if(!molecule->atoms) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + molecule->n_atoms * sizeof(struct trg_atom), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + atom = molecule->atoms; + + for(j=molecule->n_chains; j--;) + { + chain->molecule = molecule; + + memcpy(&chain->id, block->block_contents+offset, + sizeof(chain->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &chain->id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(chain->id); + + len = min(strlen(block->block_contents+offset) + 1, + TRG_MAX_STR_LEN); + chain->name = (char *) malloc(len); + strncpy(chain->name, + block->block_contents+offset, len); + offset += len; + + memcpy(&chain->n_residues, block->block_contents+offset, + sizeof(chain->n_residues)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &chain->n_residues) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(chain->n_residues); + + chain->residues = residue; + for(k=chain->n_residues; k--;) + { + residue->chain = chain; + memcpy(&residue->id, block->block_contents+offset, + sizeof(residue->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &residue->id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(residue->id); + + len = min(strlen(block->block_contents+offset) + 1, + TRG_MAX_STR_LEN); + residue->name = (char *) malloc(len); + strncpy(residue->name, + block->block_contents+offset, len); + offset += len; + + memcpy(&residue->n_atoms, block->block_contents+offset, + sizeof(residue->n_atoms)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &residue->n_atoms) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(residue->n_atoms); + + residue->atoms = atom; + for(l=residue->n_atoms; l--;) + { + atom->residue = residue; + + memcpy(&atom->id, block->block_contents+offset, + sizeof(atom->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &atom->id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(atom->id); + + len = min(strlen(block->block_contents+offset) + 1, + TRG_MAX_STR_LEN); + atom->name = (char *) malloc(len); + strncpy(atom->name, + block->block_contents+offset, len); + offset += len; + + len = min(strlen(block->block_contents+offset) + 1, + TRG_MAX_STR_LEN); + atom->atom_type = (char *) malloc(len); + strncpy(atom->atom_type, + block->block_contents+offset, len); + offset += len; + + atom++; + } + residue++; + } + chain++; + } + + memcpy(&molecule->n_bonds, block->block_contents+offset, + sizeof(molecule->n_bonds)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &molecule->n_bonds) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_bonds); + + molecule->bonds = (struct trg_bond *) malloc(molecule->n_bonds * + sizeof(struct trg_bond)); + if(!molecule->bonds) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + molecule->n_bonds * sizeof(struct trg_bond), + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + bond = molecule->bonds; + + for(j=molecule->n_bonds; j--;) + { + memcpy(&bond->from_atom_id, block->block_contents+offset, + sizeof(bond->from_atom_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &bond->from_atom_id) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(bond->from_atom_id); + + memcpy(&bond->to_atom_id, block->block_contents+offset, + sizeof(bond->to_atom_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &bond->to_atom_id) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(bond->to_atom_id); + + bond++; + } + } + + return(TRG_SUCCESS); +} + +/* FIXME: Update this according to the new specs */ +static trg_function_status trg_write_molecules_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + write_mode mode) +{ + int len = 0; + int i, j, k, l, offset = 0; + struct trg_molecule *molecule; + struct trg_chain *chain; + struct trg_residue *residue; + struct trg_atom *atom; + struct trg_bond *bond; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(mode != TRG_COPY_EXISTING && !trg_data->molecules) + { + return(TRG_SUCCESS); + } + +// printf("Predicting block size\n"); + /* First predict the size of the block */ + for(i = 0; i < trg_data->n_molecules; i++) + { + molecule = &trg_data->molecules[i]; +// printf("mol %s\n", molecule->name); + if(!molecule->name) + { + molecule->name = (char *) malloc(1); + if(!molecule->name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + molecule->name[0] = 0; + } + len += min(strlen(molecule->name) + 1, TRG_MAX_STR_LEN); + + chain = molecule->chains; + for(j = molecule->n_chains; j--;) + { + len += sizeof(chain->id); + + if(!chain->name) + { + chain->name = (char *) malloc(1); + if(!chain->name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + chain->name[0] = 0; + } + len += min(strlen(chain->name) + 1, TRG_MAX_STR_LEN); + + len += sizeof(chain->n_residues); + + chain++; + } + + residue = molecule->residues; + for(j = molecule->n_residues; j--;) + { + len += sizeof(residue->id); + + if(!residue->name) + { + residue->name = (char *) malloc(1); + if(!residue->name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + residue->name[0] = 0; + } + len += min(strlen(residue->name) + 1, TRG_MAX_STR_LEN); + + len += sizeof(residue->n_atoms); + + residue++; + } + + atom = molecule->atoms; + for(j = molecule->n_atoms; j--;) + { + len += sizeof(atom->id); + if(!atom->name) + { + atom->name = (char *) malloc(1); + if(!atom->name) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + atom->name[0] = 0; + } + len += min(strlen(atom->name) + 1, TRG_MAX_STR_LEN); + + if(!atom->atom_type) + { + atom->atom_type = (char *) malloc(1); + if(!atom->atom_type) + { + printf("Cannot allocate memory (1 byte). %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + atom->atom_type[0] = 0; + } + len += min(strlen(atom->atom_type) + 1, TRG_MAX_STR_LEN); + + atom++; + } + + bond = molecule->bonds; + for(j = molecule->n_bonds; j--;) + { + len += sizeof(bond->from_atom_id) + sizeof(bond->to_atom_id); + } + } + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + block->block_contents_size = sizeof(trg_data->n_molecules) + + (sizeof(molecule->id) + + sizeof(molecule->quaternary_str) + + sizeof(molecule->n_chains) + + sizeof(molecule->n_residues) + + sizeof(molecule->n_atoms) + + sizeof(molecule->n_bonds)) * + trg_data->n_molecules + + len; + + if(!trg_data->var_num_atoms_flag) + { + block->block_contents_size += trg_data->n_molecules * sizeof(uint64_t); + } + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + memcpy(block->block_contents+offset, &trg_data->n_molecules, + sizeof(trg_data->n_molecules)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(trg_data->n_molecules); + + for(i = 0; i < trg_data->n_molecules; i++) + { + molecule = &trg_data->molecules[i]; +// printf("i=%d\n", i); + memcpy(block->block_contents+offset, &molecule->id, + sizeof(molecule->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->id); + +// printf("Wrote id: %"PRId64" offset: %d\n", molecule->id, offset); + len = min(strlen(molecule->name) + 1, TRG_MAX_STR_LEN); + strncpy(block->block_contents + offset, molecule->name, len); + offset += len; + + memcpy(block->block_contents+offset, &molecule->quaternary_str, + sizeof(molecule->quaternary_str)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->quaternary_str); + + if(!trg_data->var_num_atoms_flag) + { + memcpy(block->block_contents+offset, + &trg_data->molecule_cnt_list[i], sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + } + + memcpy(block->block_contents+offset, &molecule->n_chains, + sizeof(molecule->n_chains)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_chains); + + memcpy(block->block_contents+offset, &molecule->n_residues, + sizeof(molecule->n_residues)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_residues); + + memcpy(block->block_contents+offset, &molecule->n_atoms, + sizeof(molecule->n_atoms)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_atoms); + + chain = molecule->chains; + for(j = molecule->n_chains; j--;) + { + memcpy(block->block_contents+offset, &chain->id, sizeof(chain->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(chain->id); + + len = min(strlen(chain->name) + 1, TRG_MAX_STR_LEN); + strncpy(block->block_contents + offset, chain->name, len); + offset += len; + + memcpy(block->block_contents+offset, &chain->n_residues, + sizeof(chain->n_residues)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(chain->n_residues); + + residue = chain->residues; + for(k = chain->n_residues; k--;) + { + memcpy(block->block_contents+offset, &residue->id, sizeof(residue->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(residue->id); + + len = min(strlen(residue->name) + 1, TRG_MAX_STR_LEN); + strncpy(block->block_contents + offset, residue->name, len); + offset += len; + + memcpy(block->block_contents+offset, &residue->n_atoms, + sizeof(residue->n_atoms)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(residue->n_atoms); + + atom = residue->atoms; + for(l = residue->n_atoms; l--;) + { + // printf("j=%d\n", j); + memcpy(block->block_contents+offset, &atom->id, sizeof(atom->id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(atom->id); + + len = min(strlen(atom->name) + 1, TRG_MAX_STR_LEN); + strncpy(block->block_contents + offset, atom->name, len); + offset += len; + + len = min(strlen(atom->atom_type) + 1, TRG_MAX_STR_LEN); + strncpy(block->block_contents + offset, atom->atom_type, len); + offset += len; + + atom++; + } + residue++; + } + chain++; + } + + memcpy(block->block_contents+offset, &molecule->n_bonds, + sizeof(molecule->n_bonds)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(molecule->n_bonds); + + bond = molecule->bonds; + for(j = molecule->n_bonds; j--;) + { + memcpy(block->block_contents+offset, &bond->from_atom_id, + sizeof(bond->from_atom_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(bond->from_atom_id); + + memcpy(block->block_contents+offset, &bond->to_atom_id, + sizeof(bond->to_atom_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(bond->to_atom_id); + + bond++; + } + } + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + + +static trg_function_status trg_read_frame_set_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block) +{ + int i, file_pos, offset = 0; + uint64_t prev_n_particles; + trg_bool same_hash; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + struct trg_particle_mapping *mapping; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to + * disk even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("Frame set block contents corrupt. Hashes do not match. " + "%s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + file_pos = ftell(trg_data->input_file); + + if(frame_set->n_mapping_blocks && frame_set->mappings) + { + for(i = frame_set->n_mapping_blocks; i--;) + { + mapping = &frame_set->mappings[i]; + if(mapping->real_particle_numbers) + { + free(mapping->real_particle_numbers); + mapping->real_particle_numbers = 0; + } + } + free(frame_set->mappings); + frame_set->n_mapping_blocks = 0; + } + + if(trg_data->first_trajectory_frame_set_input_file_pos <= 0) + { + trg_data->first_trajectory_frame_set_input_file_pos = file_pos; + } + /* FIXME: Should check the frame number instead of the file_pos, in case + * frame sets are not in order */ + if(trg_data->last_trajectory_frame_set_input_file_pos < file_pos) + { + trg_data->last_trajectory_frame_set_input_file_pos = file_pos; + } + + memcpy(&frame_set->first_frame, block->block_contents, + sizeof(frame_set->first_frame)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &frame_set->first_frame) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->first_frame); + + memcpy(&frame_set->n_frames, block->block_contents + offset, + sizeof(frame_set->n_frames)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &frame_set->n_frames) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->n_frames); + + if(trg_data->var_num_atoms_flag) + { + prev_n_particles = frame_set->n_particles; + frame_set->n_particles = 0; + /* If the list of molecule counts has already been created assume that + * it is of correct size. */ + if(!frame_set->molecule_cnt_list) + { + frame_set->molecule_cnt_list = + (uint64_t *) malloc(sizeof(uint64_t) * trg_data->n_molecules); + + if(!frame_set->molecule_cnt_list) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(uint64_t) * trg_data->n_molecules, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + } + for(i = 0; i < trg_data->n_molecules; i++) + { + memcpy(&frame_set->molecule_cnt_list[i], + block->block_contents + offset, + sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set->molecule_cnt_list[i]) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + frame_set->n_particles += trg_data->molecules[i].n_atoms * + frame_set->molecule_cnt_list[i]; + } + if(prev_n_particles && frame_set->n_particles != prev_n_particles) + { + /* FIXME: Particle dependent data memory management */ + } + } + + memcpy(&frame_set->next_frame_set_file_pos, + block->block_contents + offset, + sizeof(frame_set->next_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set->next_frame_set_file_pos) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->next_frame_set_file_pos); + + memcpy(&frame_set->prev_frame_set_file_pos, + block->block_contents + offset, + sizeof(frame_set->prev_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set->prev_frame_set_file_pos) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->prev_frame_set_file_pos); + + memcpy(&frame_set->long_stride_next_frame_set_file_pos, + block->block_contents + offset, + sizeof(frame_set->long_stride_next_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set-> + long_stride_next_frame_set_file_pos) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->long_stride_next_frame_set_file_pos); + + memcpy(&frame_set->long_stride_prev_frame_set_file_pos, + block->block_contents + offset, + sizeof(frame_set->long_stride_prev_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set-> + long_stride_prev_frame_set_file_pos) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->long_stride_prev_frame_set_file_pos); + +// frame_set->n_double_particle_data_blocks = 0; +// frame_set->n_float_particle_data_blocks = 0; +// frame_set->n_int_particle_data_blocks = 0; +// frame_set->n_gen_particle_data_blocks = 0; +// frame_set->n_double_data_blocks = 0; +// frame_set->n_float_data_blocks = 0; +// frame_set->n_int_data_blocks = 0; +// frame_set->n_gen_data_blocks = 0; + + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_frame_set_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + write_mode mode) +{ + char *temp_name; + uint64_t i; + int offset = 0, name_len; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + name_len = strlen("TRAJECTORY FRAME SET"); + + if(!block->name || strlen(block->name) < name_len) + { + temp_name = realloc(block->name, name_len + 1); + if(!temp_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + name_len+1, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + block->name = temp_name; + } + strcpy(block->name, "TRAJECTORY FRAME SET"); + + block->block_contents_size = sizeof(uint64_t) * 6; + if(trg_data->var_num_atoms_flag) + { + block->block_contents_size += sizeof(uint64_t) * trg_data->n_molecules; + } + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + memcpy(block->block_contents, &frame_set->first_frame, + sizeof(frame_set->first_frame)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->first_frame); + + memcpy(block->block_contents+offset, &frame_set->n_frames, + sizeof(frame_set->n_frames)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->n_frames); + + if(trg_data->var_num_atoms_flag) + { + for(i = 0; i < trg_data->n_molecules; i++) + { + memcpy(block->block_contents+offset, + &frame_set->molecule_cnt_list[i], + sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + } + } + + + memcpy(block->block_contents+offset, &frame_set->next_frame_set_file_pos, + sizeof(frame_set->next_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->next_frame_set_file_pos); + + memcpy(block->block_contents+offset, &frame_set->prev_frame_set_file_pos, + sizeof(frame_set->prev_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->prev_frame_set_file_pos); + + memcpy(block->block_contents+offset, + &frame_set->long_stride_next_frame_set_file_pos, + sizeof(frame_set->long_stride_next_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->long_stride_next_frame_set_file_pos); + + memcpy(block->block_contents+offset, + &frame_set->long_stride_prev_frame_set_file_pos, + sizeof(frame_set->long_stride_prev_frame_set_file_pos)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(frame_set->long_stride_prev_frame_set_file_pos); + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_trajectory_toc_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block) +{ + uint64_t i, old_n_blocks; + int offset = 0, len; + trg_bool same_hash; + struct trg_frame_set_toc *toc = + &trg_data->current_trajectory_frame_set.contents; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to disk + * even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("Table of contents block contents corrupt. Hashes do not match. " + "%s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + old_n_blocks = toc->n_blocks; + + memcpy(&toc->n_blocks, block->block_contents, + sizeof(toc->n_blocks)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &toc->n_blocks) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(toc->n_blocks); + + if(old_n_blocks != toc->n_blocks) + { + if(toc->block_names) + { + for(i = old_n_blocks; i--;) + { + if(toc->block_names[i]) + { + free(toc->block_names[i]); + } + } + free(toc->block_names); + } + toc->block_names = (char **) malloc(toc->n_blocks * sizeof(char *)); + if(!toc->block_names) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + toc->n_blocks * sizeof(uint64_t), __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + } + + for(i = 0; i < toc->n_blocks; i++) + { + len = min(strlen(block->block_contents+offset) + 1, TRG_MAX_STR_LEN); + + toc->block_names[i] = (char *) malloc(len); + if(!toc->block_names[i]) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + len, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + strncpy(toc->block_names[i], block->block_contents+offset, len); + + offset += len; + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_trajectory_toc_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + write_mode mode) +{ + char *temp_name; + uint64_t i; + int offset = 0, name_len; + struct trg_frame_set_toc *toc = + &trg_data->current_trajectory_frame_set.contents; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + name_len = strlen("BLOCK TABLE OF CONTENTS"); + + if(!block->name || strlen(block->name) < name_len) + { + temp_name = realloc(block->name, name_len + 1); + if(!temp_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + name_len+1, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + block->name = temp_name; + } + strcpy(block->name, "BLOCK TABLE OF CONTENTS"); + + block->block_contents_size = sizeof(uint64_t); + + for(i = toc->n_blocks; i--;) + { + block->block_contents_size += strlen(toc->block_names[i]) + 1; + } + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + memcpy(block->block_contents, &toc->n_blocks, sizeof(toc->n_blocks)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(toc->n_blocks); + + for(i = 0; i < toc->n_blocks; i++) + { + name_len = strlen(toc->block_names[i]) + 1; + strncpy(block->block_contents+offset, toc->block_names[i], name_len); + offset += name_len; + } + + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_trajectory_mapping_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block) +{ + uint64_t i, old_n_particles; + int offset = 0; + trg_bool same_hash; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_particle_mapping *mapping, *mappings; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to disk + * even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("Particle mapping block contents corrupt. Hashes do not match. " + "%s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + frame_set->n_mapping_blocks++; + mappings = realloc(frame_set->mappings, + sizeof(struct trg_particle_mapping) * + frame_set->n_mapping_blocks); + if(!mappings) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + frame_set->mappings = mappings; + mapping = &mappings[frame_set->n_mapping_blocks - 1]; + + + memcpy(&mapping->num_first_particle, block->block_contents+offset, + sizeof(mapping->num_first_particle)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &mapping->num_first_particle) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(mapping->num_first_particle); + + old_n_particles = mapping->n_particles; + + memcpy(&mapping->n_particles, block->block_contents+offset, + sizeof(mapping->n_particles)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &mapping->n_particles) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(mapping->n_particles); + + if(old_n_particles != mapping->n_particles) + { + if(mapping->real_particle_numbers) + { + free(mapping->real_particle_numbers); + } + mapping->real_particle_numbers = malloc(mapping->n_particles * + sizeof(uint64_t)); + if(!mapping->real_particle_numbers) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + mapping->n_particles * sizeof(uint64_t), __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + } + + for(i = 0; i < mapping->n_particles; i++) + { + memcpy(&mapping->real_particle_numbers[i], + block->block_contents + offset, + sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &mapping->real_particle_numbers[i]) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_trajectory_mapping_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + int mapping_block_nr, + write_mode mode) +{ + char *temp_name; + int i, offset = 0, name_len; + struct trg_particle_mapping *mapping = + &trg_data->current_trajectory_frame_set.mappings[mapping_block_nr]; + + if(mapping_block_nr >= + trg_data->current_trajectory_frame_set.n_mapping_blocks) + { + printf("Mapping block index out of bounds. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + name_len = strlen("PARTICLE MAPPING"); + + if(!block->name || strlen(block->name) < name_len) + { + temp_name = realloc(block->name, name_len + 1); + if(!temp_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + name_len+1, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + block->name = temp_name; + } + strcpy(block->name, "PARTICLE MAPPING"); + + block->block_contents_size = sizeof(uint64_t) * (2 + mapping->n_particles); + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + memcpy(block->block_contents, &mapping->num_first_particle, + sizeof(mapping->num_first_particle)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(mapping->num_first_particle); + + memcpy(block->block_contents+offset, &mapping->n_particles, + sizeof(mapping->n_particles)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(mapping->n_particles); + + for(i = 0; i < mapping->n_particles; i++) + { + memcpy(block->block_contents+offset, &mapping->real_particle_numbers[i], + sizeof(uint64_t)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(uint64_t); + } + + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_create_particle_data_block + (struct trg_trajectory *trg_data, + const trg_block_type block_type_flag) +{ + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_particle_data *data; + + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + frame_set->n_particle_data_blocks++; + data = realloc(frame_set->tr_particle_data, + sizeof(struct trg_particle_data) * + frame_set->n_particle_data_blocks); + if(!data) + { + printf("Cannot allocate memory (%lu bytes). %s: %d\n", + sizeof(struct trg_particle_data) * + frame_set->n_particle_data_blocks, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + frame_set->tr_particle_data = data; + data = &frame_set->tr_particle_data[frame_set-> + n_particle_data_blocks - 1]; + } + else + { + trg_data->n_particle_data_blocks++; + data = realloc(trg_data->non_tr_particle_data, + sizeof(struct trg_particle_data) * + trg_data->n_particle_data_blocks); + if(!data) + { + printf("Cannot allocate memory (%lu bytes). %s: %d\n", + sizeof(struct trg_particle_data) * + trg_data->n_particle_data_blocks, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + trg_data->non_tr_particle_data = data; + data = &trg_data->non_tr_particle_data[trg_data-> + n_particle_data_blocks - 1]; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_allocate_particle_data_mem + (struct trg_trajectory *trg_data, + struct trg_particle_data *data, + uint64_t n_frames, + const uint64_t n_particles, + const uint64_t n_values_per_frame) +{ + union data_values ***values; + uint64_t i, j, k; + + for(i = data->n_frames; i--;) + { + for(j = n_particles; j--;) + { + free(data->values[i][j]); + } + free(data->values[i]); + } + data->n_frames = n_frames; + n_frames = max(1, n_frames); + data->n_values_per_frame = n_values_per_frame; + values = (union data_values ***) realloc(data->values, + sizeof(union data_values **) * + n_frames); + if(!values) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(union data_values **) * n_frames, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + data->values = values; + + for(i = n_frames; i-- ;) + { + data->values[i] = (union data_values **) + malloc(sizeof(union data_values *) * + n_particles); + if(!data->values[i]) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(union data_values *) * n_particles, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + for(j = n_particles; j--;) + { + data->values[i][j] = (union data_values *) + malloc(sizeof(union data_values) * + n_values_per_frame); + if(!data->values[i][j]) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(union data_values) * n_values_per_frame, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(data->datatype == TRG_CHAR_DATA) + { + for(k = n_values_per_frame; k--;) + { + data->values[i][j][k].c = 0; + } + } + } + } + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_particle_data + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + int *offset, + const char datatype, + const uint64_t first_particle_number, + const uint64_t n_particles, + const uint64_t first_frame_with_data, + const uint64_t stride_length, + uint64_t n_frames, + const uint64_t n_values, + const uint64_t codec_id, + const uint64_t multiplier) +{ + uint64_t block_index, i, j, k, tot_n_particles; + int size, len; + struct trg_particle_data *data; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + trg_block_type block_type_flag; + + switch(datatype) + { + case TRG_CHAR_DATA: + size = 1; + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + } + + if(trg_data->current_trajectory_frame_set_input_file_pos > 0) + { + block_type_flag = TRG_TRAJECTORY_BLOCK; + } + else + { + block_type_flag = TRG_NON_TRAJECTORY_BLOCK; + } + + block_index = -1; + /* See if there is already a data block of this ID */ + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + for(i = frame_set->n_particle_data_blocks; i-- ;) + { + data = &frame_set->tr_particle_data[i]; + if(data->block_id == block->id) + { + block_index = i; + break; + } + } + } + else + { + for(i = trg_data->n_particle_data_blocks; i-- ;) + { + data = &trg_data->non_tr_particle_data[i]; + if(data->block_id == block->id) + { + block_index = i; + break; + } + } + } + + /* Otherwise create a data block */ + if(block_index == -1) + { + if(trg_create_particle_data_block(trg_data, block_type_flag) != + TRG_SUCCESS) + { + printf("Cannot create particle data block. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + data = &frame_set->tr_particle_data[frame_set-> + n_particle_data_blocks - 1]; + } + else + { + data = &trg_data->non_tr_particle_data[trg_data-> + n_particle_data_blocks - 1]; + } + data->block_id = block->id; + + data->block_name = malloc(strlen(block->name) + 1); + if(!data->block_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + (int)strlen(block->name)+1, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strcpy(data->block_name, block->name); + + data->datatype = datatype; + + data->values = 0; + data->n_frames = 0; + data->codec_id = codec_id; + data->compression_multiplier = multiplier; + } + + if(block_type_flag == TRG_TRAJECTORY_BLOCK && + trg_data->var_num_atoms_flag) + { + tot_n_particles = frame_set->n_particles; + } + else + { + tot_n_particles = trg_data->n_particles; + } + + /* Allocate memory */ + if(!data->values || data->n_frames != n_frames || + data->n_values_per_frame != n_values) + { + if(trg_allocate_particle_data_mem(trg_data, data, n_frames, + tot_n_particles, n_values) != + TRG_SUCCESS) + { + printf("Cannot allocate memory for particle data. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + data->first_frame_with_data = first_frame_with_data; + data->stride_length = stride_length; + + for(i = 0; i < n_frames; i++) + { + for(j = first_particle_number; j < n_particles; j++) + { + for(k = 0; k < n_values; k++) + { + /* FIXME: Not optimal to have this switch statement inside + * the for loops. Check if it should be moved outside. */ + switch(datatype) + { + case TRG_FLOAT_DATA: + memcpy(&data->values[i][j][k].f, + block->block_contents+*offset, + size); + if(trg_data->endianness_32 != TRG_BIG_ENDIAN_32) + { + if(trg_swap_byte_order_32(trg_data, + (uint32_t *) &data->values[i][j][k]) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + break; + case TRG_INT_DATA: + memcpy(&data->values[i][j][k].i, + block->block_contents+*offset, + size); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *) &data->values[i][j][k].i) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + break; + case TRG_CHAR_DATA: + len = min(strlen(block->block_contents+*offset) + 1, + TRG_MAX_STR_LEN); + if(data->values[i][j][k].c) + { + free(data->values[i][j][k].c); + } + data->values[i][j][k].c = malloc(len); + if(!data->values[i][j][k].c) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + len, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->values[i][j][k].c, + block->block_contents+*offset, len); + *offset += len; + break; + case TRG_DOUBLE_DATA: + default: + memcpy(&data->values[i][j][k].d, + block->block_contents+*offset, + size); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *) &data->values[i][j][k].d) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + } + } + } + } + return(TRG_SUCCESS); +} + + +static trg_function_status trg_write_particle_data_block + (struct trg_trajectory *trg_data, + struct trg_gen_block *block, + const int block_index, + const struct trg_particle_mapping *mapping, + const write_mode mode) +{ + uint64_t n_particles, num_first_particle, n_frames; + int i, j, k, offset = 0, size, len; + char temp, *temp_name; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_particle_data *data; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(trg_data->current_trajectory_frame_set_output_file_pos > 0) + { + data = &frame_set->tr_particle_data[block_index]; + } + else + { + data = &trg_data->non_tr_particle_data[block_index]; + } + + switch(data->datatype) + { + case TRG_CHAR_DATA: + size = 1; + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + } + + len = strlen(data->block_name) + 1; + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + if(!block->name || strlen(block->name) < len) + { + temp_name = realloc(block->name, len); + if(!temp_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + block->name = temp_name; + } + strncpy(block->name, data->block_name, len); + + /* If writing frame independent data data->n_frames is be 0, but n_frames + is used for the loop writing the data (and reserving memory) and needs + to be at least 1 */ + n_frames = max(1, data->n_frames); + + + if(mapping && mapping->n_particles != 0) + { + n_particles = mapping->n_particles; + num_first_particle = mapping->num_first_particle; + } + else + { + num_first_particle = 0; + if(trg_data->var_num_atoms_flag) + { + n_particles = frame_set->n_particles; + } + else + { + n_particles = trg_data->n_particles; + } + } + + block->block_contents_size = sizeof(char) * 3 + + sizeof(data->n_values_per_frame) + + sizeof(data->codec_id) + + sizeof(num_first_particle) + + sizeof(n_particles); + + if(data->codec_id != TRG_UNCOMPRESSED) + { + block->block_contents_size += sizeof(data->compression_multiplier); + } + + if(data->datatype == TRG_CHAR_DATA) + { + for(i = n_frames; i--;) + { + for(j = num_first_particle; j < n_particles; j++) + { + for(k = data->n_values_per_frame; k--;) + { + block->block_contents_size += + strlen(data->values[i][j][k].c) + 1; + } + } + } + } + else + { + block->block_contents_size += size * n_frames * n_particles * + data->n_values_per_frame; + } + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + + memcpy(block->block_contents, &data->datatype, sizeof(char)); + offset += sizeof(char); + + if(data->n_frames > 0) + { + temp = TRG_FRAME_DEPENDENT + TRG_PARTICLE_DEPENDENT; + } + else + { + temp = TRG_PARTICLE_DEPENDENT; + } + memcpy(block->block_contents+offset, &temp, sizeof(char)); + offset += sizeof(char); + + if(data->n_frames > 0 && data->stride_length > 1) + { + temp = 1; + } + else + { + temp = 0; + } + memcpy(block->block_contents+offset, &temp, sizeof(char)); + offset += sizeof(char); + + memcpy(block->block_contents+offset, &data->n_values_per_frame, + sizeof(data->n_values_per_frame)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->n_values_per_frame); + + memcpy(block->block_contents+offset, &data->codec_id, + sizeof(data->codec_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->codec_id); + + if(data->codec_id != TRG_UNCOMPRESSED) + { + memcpy(block->block_contents+offset, &data->compression_multiplier, + sizeof(data->compression_multiplier)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->compression_multiplier); + } + + if(data->n_frames > 0 && data->stride_length > 1) + { + memcpy(block->block_contents+offset, &data->first_frame_with_data, + sizeof(data->first_frame_with_data)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->first_frame_with_data); + + memcpy(block->block_contents+offset, &data->stride_length, + sizeof(data->stride_length)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->stride_length); + } + + + memcpy(block->block_contents+offset, &num_first_particle, + sizeof(num_first_particle)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(num_first_particle); + + memcpy(block->block_contents+offset, &n_particles, sizeof(n_particles)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents+offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(n_particles); + + + for(i = 0; i < data->n_frames; i++) + { + for(j = num_first_particle; j < n_particles; j++) + { + for(k = 0; k < data->n_values_per_frame; k++) + { + switch(data->datatype) + { + case TRG_FLOAT_DATA: + memcpy(block->block_contents+offset, + &data->values[i][j][k].f, + size); + break; + case TRG_INT_DATA: + memcpy(block->block_contents+offset, + &data->values[i][j][k].i, + size); + break; + case TRG_CHAR_DATA: + len = strlen(data->values[i][j][k].c) + 1; + strncpy(block->block_contents+offset, + data->values[i][j][k].c, len); + offset += len; + break; + case TRG_DOUBLE_DATA: + default: + memcpy(block->block_contents+offset, + &data->values[i][j][k].d, + size); + } + if(data->datatype == TRG_FLOAT_DATA) + { + if(trg_data->endianness_32 != TRG_BIG_ENDIAN_32) + { + if(trg_swap_byte_order_32(trg_data, + (uint32_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(size); + } + else if(data->datatype == TRG_INT_DATA || + data->datatype == TRG_DOUBLE_DATA) + { + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(size); + } + } + } + } + + + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +/* UNTESTED */ +static trg_function_status trg_create_data_block + (struct trg_trajectory *trg_data, + const trg_block_type block_type_flag) +{ + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_data *data; + + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + frame_set->n_data_blocks++; + data = realloc(frame_set->tr_data, sizeof(struct trg_data) * + frame_set->n_data_blocks); + if(!data) + { + printf("Cannot allocate memory (%lu bytes). %s: %d\n", + sizeof(struct trg_data) * frame_set->n_data_blocks, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + frame_set->tr_data = data; + data = &frame_set->tr_data[frame_set->n_data_blocks - 1]; + } + else + { + trg_data->n_data_blocks++; + data = realloc(trg_data->non_tr_data, sizeof(struct trg_data) * + trg_data->n_data_blocks); + if(!data) + { + printf("Cannot allocate memory (%lu bytes). %s: %d\n", + sizeof(struct trg_data) * trg_data->n_data_blocks, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + trg_data->non_tr_data = data; + data = &trg_data->non_tr_data[trg_data->n_data_blocks - 1]; + } + + return(TRG_SUCCESS); +} + +/* UNTESTED */ +trg_function_status trg_allocate_data_mem + (struct trg_trajectory *trg_data, + struct trg_data *data, + uint64_t n_frames, + const uint64_t n_values_per_frame) +{ + union data_values **values; + uint64_t i, j; + + for(i = data->n_frames; i--;) + { + if(data->datatype == TRG_CHAR_DATA) + { + for(j = data->n_values_per_frame; j--;) + { + if(data->values[i][j].c) + { + free(data->values[i][j].c); + data->values[i][j].c = 0; + } + } + } + free(data->values[i]); + } + data->n_frames = n_frames; + n_frames = max(1, n_frames); + data->n_values_per_frame = n_values_per_frame; + values = (union data_values **) realloc(data->values, + sizeof(union data_values *) * + n_frames); + if(!values) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(union data_values **) * n_frames, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + data->values = values; + + for(i = n_frames; i-- ;) + { + data->values[i] = (union data_values *) + malloc(sizeof(union data_values) * + n_values_per_frame); + if(!data->values[i]) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(union data_values) * n_values_per_frame, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(data->datatype == TRG_CHAR_DATA) + { + for(j = n_values_per_frame; j--;) + { + data->values[i][j].c = 0; + } + } + } + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_data(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + int *offset, + const char datatype, + const uint64_t first_frame_with_data, + const uint64_t stride_length, + uint64_t n_frames, + const uint64_t n_values, + const uint64_t codec_id, + const uint64_t multiplier) +{ + uint64_t block_index, i, j; + int size, len; + struct trg_data *data; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + trg_block_type block_type_flag; + +// printf("%s\n", block->name); + + if(trg_data->current_trajectory_frame_set_input_file_pos > 0) + { + block_type_flag = TRG_TRAJECTORY_BLOCK; + } + else + { + block_type_flag = TRG_NON_TRAJECTORY_BLOCK; + } + + switch(datatype) + { + case TRG_CHAR_DATA: + size = 1; + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + } + + block_index = -1; + /* See if there is already a data block of this ID */ + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + for(i = frame_set->n_data_blocks; i-- ;) + { + data = &frame_set->tr_data[i]; + if(data->block_id == block->id) + { + block_index = i; + break; + } + } + } + else + { + for(i = trg_data->n_data_blocks; i-- ;) + { + data = &trg_data->non_tr_data[i]; + if(data->block_id == block->id) + { + block_index = i; + break; + } + } + } + + /* Otherwise create a data block */ + if(block_index == -1) + { + if(trg_create_data_block(trg_data, block_type_flag) != + TRG_SUCCESS) + { + printf("Cannot create particle data block. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + data = &frame_set->tr_data[frame_set->n_data_blocks - 1]; + } + else + { + data = &trg_data->non_tr_data[trg_data->n_data_blocks - 1]; + } + data->block_id = block->id; + + data->block_name = malloc(strlen(block->name) + 1); + if(!data->block_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + (int)strlen(block->name)+1, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + strcpy(data->block_name, block->name); + + data->datatype = datatype; + + data->values = 0; + data->n_frames = 0; + data->codec_id = codec_id; + data->compression_multiplier = multiplier; + } + + /* Allocate memory */ + if(!data->values || data->n_frames != n_frames || + data->n_values_per_frame != n_values) + { + if(trg_allocate_data_mem(trg_data, data, n_frames, n_values) != + TRG_SUCCESS) + { + printf("Cannot allocate memory for data. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + data->first_frame_with_data = first_frame_with_data; + data->stride_length = stride_length; + + n_frames = max(1, n_frames); + + for(i = 0; i < n_frames; i++) + { + for(j = 0; j < n_values; j++) + { + switch(datatype) + { + case TRG_CHAR_DATA: + len = min(strlen(block->block_contents+*offset) + 1, + TRG_MAX_STR_LEN); + if(data->values[i][j].c) + { + free(data->values[i][j].c); + } + data->values[i][j].c = malloc(len); + if(!data->values[i][j].c) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + len, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->values[i][j].c, block->block_contents+*offset, + len); + *offset += len; + break; + case TRG_INT_DATA: + memcpy(&data->values[i][j].i, block->block_contents+*offset, + size); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *) &data->values[i][j].i) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + break; + case TRG_FLOAT_DATA: + memcpy(&data->values[i][j].f, block->block_contents+*offset, + size); + if(trg_data->endianness_32 != TRG_BIG_ENDIAN_32) + { + if(trg_swap_byte_order_32(trg_data, + (uint32_t *) &data->values[i][j]) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + break; + case TRG_DOUBLE_DATA: + default: + memcpy(&data->values[i][j].d, block->block_contents+*offset, + size); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *) &data->values[i][j].d) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + *offset += size; + } + } + } + return(TRG_SUCCESS); +} + +static trg_function_status trg_write_data_block(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + const int block_index, + const write_mode mode) +{ + uint64_t n_frames; + int i, j, offset = 0, size, len; + char temp, *temp_name; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_data *data; + trg_block_type block_type_flag; + + if(trg_data->current_trajectory_frame_set_output_file_pos > 0) + { + block_type_flag = TRG_TRAJECTORY_BLOCK; + } + else + { + block_type_flag = TRG_NON_TRAJECTORY_BLOCK; + } + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + data = &frame_set->tr_data[block_index]; + } + else + { + data = &trg_data->non_tr_data[block_index]; + } + + switch(data->datatype) + { + case TRG_CHAR_DATA: + size = 1; + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + } + + /* If just dumping the whole block_contents it is not certain that the + * contents are known beforehand (e.g. due to different file versions) */ + if(mode == TRG_COPY_EXISTING) + { + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(!block->block_contents) + { + printf("No block data to write. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + return(TRG_SUCCESS); + } + + len = strlen(data->block_name) + 1; + + if(!block->name || strlen(block->name) < len) + { + temp_name = realloc(block->name, len); + if(!temp_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len+1, + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + block->name = temp_name; + } + strncpy(block->name, data->block_name, len); + + /* If writing frame independent data data->n_frames is be 0, but n_frames + is used for the loop writing the data (and reserving memory) and needs + to be at least 1 */ + n_frames = max(1, data->n_frames); + + + block->block_contents_size = sizeof(char) * 2 + + sizeof(data->n_values_per_frame) + + sizeof(data->codec_id); + + if(data->codec_id != TRG_UNCOMPRESSED) + { + block->block_contents_size += sizeof(data->compression_multiplier); + } + + if(data->datatype == TRG_CHAR_DATA) + { + for(i = n_frames; i--;) + { + for(j = data->n_values_per_frame; j--;) + { + block->block_contents_size += strlen(data->values[i][j].c) + 1; + } + } + } + else + { + block->block_contents_size += size * n_frames * + data->n_values_per_frame; + } + + if(data->n_frames > 0 || data->stride_length != 0) + { + block->block_contents_size += sizeof(char); + } + + if(block->block_contents) + { + free(block->block_contents); + } + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + + memcpy(block->block_contents, &data->datatype, sizeof(char)); + offset += sizeof(char); + + if(data->n_frames > 0 || data->stride_length != 0) + { + temp = TRG_FRAME_DEPENDENT; + } + else + { + temp = 0; + } + memcpy(block->block_contents+offset, &temp, sizeof(char)); + offset += sizeof(char); + + if(data->n_frames > 0) + { + if(data->stride_length > 1) + { + temp = 1; + } + else + { + temp = 0; + } + memcpy(block->block_contents+offset, &temp, sizeof(char)); + offset += sizeof(char); + } + + memcpy(block->block_contents+offset, &data->n_values_per_frame, + sizeof(data->n_values_per_frame)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->n_values_per_frame); + + memcpy(block->block_contents+offset, &data->codec_id, + sizeof(data->codec_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->codec_id); + + if(data->codec_id != TRG_UNCOMPRESSED) + { + memcpy(block->block_contents+offset, &data->compression_multiplier, + sizeof(data->compression_multiplier)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->compression_multiplier); + } + + if(data->n_frames > 0 && data->stride_length > 1) + { + memcpy(block->block_contents+offset, &data->first_frame_with_data, + sizeof(data->first_frame_with_data)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->first_frame_with_data); + + memcpy(block->block_contents+offset, &data->stride_length, + sizeof(data->stride_length)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *) + (block->block_contents + offset)) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(data->stride_length); + } + + for(i = 0; i < n_frames; i++) + { + for(j = 0; j < data->n_values_per_frame; j++) + { + switch(data->datatype) + { + case TRG_CHAR_DATA: + len = strlen(data->values[i][j].c) + 1; + strncpy(block->block_contents+offset, data->values[i][j].c, + len); + offset += len; + break; + case TRG_INT_DATA: + memcpy(block->block_contents+offset, &data->values[i][j].i, + size); + break; + case TRG_FLOAT_DATA: + memcpy(block->block_contents+offset, &data->values[i][j].f, + size); + break; + case TRG_DOUBLE_DATA: + default: + memcpy(block->block_contents+offset, &data->values[i][j].d, + size); + } + if(data->datatype != TRG_CHAR_DATA) + { + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + (uint64_t *)(block->block_contents+offset)) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. " + "%s: %d\n", + __FILE__, __LINE__); + } + } + offset += size; + } + } + } + + if(trg_write_block_header(trg_data, block, mode) != TRG_SUCCESS) + { + printf("Cannot write header of file %s. %s: %d\n", + trg_data->output_file_path, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(fwrite(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) != 1) + { + printf("Could not write all block data. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_read_data_block_contents + (struct trg_trajectory *trg_data, + struct trg_gen_block *block) +{ + uint64_t n_values, codec_id, n_frames, first_frame_with_data; + uint64_t steps_between_data, block_n_particles, first_particle_number; + double multiplier; + char datatype, dependency, sparse_data; + int offset = 0; + trg_bool same_hash; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = (char *) malloc(block->block_contents_size); + if(!block->block_contents) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + block->block_contents_size, __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* Read the whole block into block_contents to be able to write it to + * disk even if it cannot be interpreted. */ + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->input_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(block); + return(TRG_CRITICAL); + } + + /* FIXME: Does not check if the size of the contents matches the expected + * size or if the contents can be read. */ + + + if(trg_verify_hash_match(block, &same_hash) != TRG_SUCCESS) + { + printf("Error comparing hashes. %s: %d\n", __FILE__, __LINE__); + return(TRG_FAILURE); + } + if(same_hash != TRUE) + { + printf("%s\n", block->hash); + printf("Data block contents corrupt. Hashes do not match. %s: %d\n", + __FILE__, __LINE__); +// return(TRG_FAILURE); + } + + memcpy(&datatype, block->block_contents+offset, + sizeof(datatype)); + offset += sizeof(datatype); + + memcpy(&dependency, block->block_contents+offset, + sizeof(dependency)); + offset += sizeof(dependency); + +// memcpy(&var_n_particles, block->block_contents+offset, +// sizeof(var_n_particles)); +// offset += sizeof(var_n_particles); + + if(dependency & TRG_FRAME_DEPENDENT) + { + memcpy(&sparse_data, block->block_contents+offset, + sizeof(sparse_data)); + offset += sizeof(sparse_data); + } + +// memcpy(&var_n_values, block->block_contents+offset, +// sizeof(var_n_values)); +// offset += sizeof(var_n_values); +// +// if(!var_n_values) +// { + memcpy(&n_values, block->block_contents+offset, + sizeof(n_values)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &n_values) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(n_values); +// } + + memcpy(&codec_id, block->block_contents+offset, + sizeof(codec_id)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &codec_id) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(codec_id); + + if(codec_id != TRG_UNCOMPRESSED) + { + memcpy(&multiplier, block->block_contents+offset, + sizeof(multiplier)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, (uint64_t *)&multiplier) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(multiplier); + } + else + { + multiplier = 1; + } + + if(dependency & TRG_FRAME_DEPENDENT) + { + if(sparse_data) + { + memcpy(&first_frame_with_data, block->block_contents+offset, + sizeof(first_frame_with_data)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &first_frame_with_data) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(first_frame_with_data); + + memcpy(&steps_between_data, block->block_contents+offset, + sizeof(steps_between_data)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &steps_between_data) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(steps_between_data); + } + else + { + first_frame_with_data = 0; + steps_between_data = 0; + } + n_frames = trg_data->current_trajectory_frame_set.n_frames; + } + else + { + first_frame_with_data = 0; + steps_between_data = 0; + n_frames = 0; + } + + if (dependency & TRG_PARTICLE_DEPENDENT) + { + memcpy(&first_particle_number, block->block_contents+offset, + sizeof(first_particle_number)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &first_particle_number) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(first_particle_number); + + memcpy(&block_n_particles, block->block_contents+offset, + sizeof(block_n_particles)); + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &block_n_particles) != + TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + offset += sizeof(block_n_particles); + } + + if (dependency & TRG_PARTICLE_DEPENDENT) + { + return(trg_read_particle_data(trg_data, block, + &offset, datatype, + first_particle_number, + block_n_particles, + first_frame_with_data, + steps_between_data, + n_frames, n_values, + codec_id, multiplier)); + } + else + { + return(trg_read_data(trg_data, block, + &offset, datatype, + first_frame_with_data, + steps_between_data, + n_frames, n_values, + codec_id, multiplier)); + } +} + +static trg_function_status trg_update_md5_hash(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + uint64_t header_start_pos, + uint64_t contents_start_pos) +{ + if(block->block_contents) + { + free(block->block_contents); + } + + block->block_contents = malloc(block->block_contents_size); + fseek(trg_data->output_file, contents_start_pos, SEEK_SET); + if(fread(block->block_contents, block->block_contents_size, 1, + trg_data->output_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + trg_generate_block_hash(block); + + fseek(trg_data->output_file, header_start_pos + 3 * sizeof(uint64_t), + SEEK_SET); + fwrite(block->hash, TRG_HASH_LEN, 1, trg_data->output_file); + + return(TRG_SUCCESS); +} + + +static trg_function_status trg_update_header_pointers + (struct trg_trajectory *trg_data) +{ + struct trg_gen_block block; + FILE *temp = trg_data->input_file; + uint64_t pos, contents_start_pos; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + printf("Cannot initialise destination file. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + trg_data->input_file = trg_data->output_file; + + trg_init_block(&block); + + trg_data->output_file_pos = ftell(trg_data->output_file); + fseek(trg_data->output_file, 0, SEEK_SET); + + if(trg_read_block_header(trg_data, &block) != TRG_SUCCESS) + { + printf("Cannot read general info header. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + contents_start_pos = ftell(trg_data->output_file); + + fseek(trg_data->output_file, block.block_contents_size - 3 * + sizeof(uint64_t), SEEK_CUR); + + trg_data->input_file = temp; + +// printf("Updating header\n"); +// printf("%ld: First frame set %ld\n", ftell(trg_data->output_file), +// trg_data->first_trajectory_frame_set_output_file_pos); + + pos = trg_data->first_trajectory_frame_set_output_file_pos; + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &pos) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + if(fwrite(&pos, sizeof(uint64_t), 1, trg_data->output_file) != 1) + { + trg_destroy_block(&block); + return(TRG_CRITICAL); + } + +// printf("%ld: Last frame set %ld\n", ftell(trg_data->output_file), +// trg_data->last_trajectory_frame_set_output_file_pos); + + pos = trg_data->last_trajectory_frame_set_output_file_pos; + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &pos) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + if(fwrite(&pos, + sizeof(uint64_t), 1, trg_data->output_file) != 1) + { + trg_destroy_block(&block); + return(TRG_CRITICAL); + } + + trg_update_md5_hash(trg_data, &block, 0, contents_start_pos); + + fseek(trg_data->output_file, trg_data->output_file_pos, SEEK_SET); + + trg_destroy_block(&block); + + return(TRG_SUCCESS); +} + +static trg_function_status trg_update_frame_set_pointers + (struct trg_trajectory *trg_data) +{ + struct trg_gen_block block; + struct trg_trajectory_frame_set *frame_set; + FILE *temp = trg_data->input_file; + uint64_t pos, header_start_pos, contents_start_pos; + + if(trg_init_output_file(trg_data, FALSE) != TRG_SUCCESS) + { + printf("Cannot initialise destination file. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + trg_init_block(&block); + trg_data->output_file_pos = ftell(trg_data->output_file); + + trg_data->input_file = trg_data->output_file; + + frame_set = &trg_data->current_trajectory_frame_set; + if(frame_set->prev_frame_set_file_pos != -1ULL && + frame_set->prev_frame_set_file_pos != 0) + { + fseek(trg_data->output_file, frame_set->prev_frame_set_file_pos, + SEEK_SET); + + header_start_pos = frame_set->prev_frame_set_file_pos; + + if(trg_read_block_header(trg_data, &block) != TRG_SUCCESS) + { + printf("Cannot read frame header. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + contents_start_pos = ftell(trg_data->output_file); + + fseek(trg_data->output_file, block.block_contents_size - 4 * + sizeof(uint64_t), SEEK_CUR); + +// printf("Updating frame set\n"); +// printf("%ld: Next frame set %ld\n", ftell(trg_data->output_file), +// trg_data->current_trajectory_frame_set_output_file_pos); + + pos = trg_data->current_trajectory_frame_set_output_file_pos; + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &pos) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + if(fwrite(&pos, sizeof(uint64_t), 1, trg_data->output_file) != 1) + { + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + trg_update_md5_hash(trg_data, &block, header_start_pos, + contents_start_pos); + + fseek(trg_data->output_file, trg_data->output_file_pos, SEEK_SET); + } + + if(frame_set->long_stride_prev_frame_set_file_pos != -1ULL && + frame_set->long_stride_prev_frame_set_file_pos != 0) + { + fseek(trg_data->output_file, + frame_set->long_stride_prev_frame_set_file_pos, + SEEK_SET); + + if(trg_read_block_header(trg_data, &block) != TRG_SUCCESS) + { + printf("Cannot read frame header. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + contents_start_pos = ftell(trg_data->output_file); + + fseek(trg_data->output_file, block.block_contents_size - 2 * + sizeof(uint64_t), SEEK_CUR); + + pos = trg_data->current_trajectory_frame_set_output_file_pos; + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, &pos) != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + if(fwrite(&pos, sizeof(uint64_t), 1, trg_data->output_file) != 1) + { + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + trg_update_md5_hash(trg_data, &block, + frame_set->long_stride_prev_frame_set_file_pos, + contents_start_pos); + } + + fseek(trg_data->output_file, trg_data->output_file_pos, SEEK_SET); + + trg_data->input_file = temp; + trg_destroy_block(&block); + + return(TRG_SUCCESS); +} + + +trg_function_status trg_set_block_name(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(block->name && strlen(block->name) < len) + { + free(block->name); + block->name = 0; + } + if(!block->name) + { + block->name = (char *) malloc(len); + if(!block->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(block->name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_init_block(struct trg_gen_block *block) +{ +// printf("In trg_init_block\n"); + + block->id = -1; +/* block->hash_type = TRG_NO_HASH; + block->hash_name = 0;*/ + memcpy(block->hash, "\0\0\0\0\0\0\0\0\0\0\0\0\0\0\0\0", TRG_HASH_LEN); + block->name = 0; + block->block_version = TRG_VERSION; + block->header_contents = 0; + block->header_contents_size = 0; + block->block_contents = 0; + block->block_contents_size = 0; + + return(TRG_SUCCESS); +} + + +trg_function_status trg_destroy_block(struct trg_gen_block *block) +{ +// printf("Destroying block\n"); +/* if(block->hash_name) + { + free(block->hash_name); + block->hash_name = 0; + }*/ +/* if(block->hash) + { + free(block->hash); + block->hash = 0; + }*/ + if(block->name) + { + free(block->name); + block->name = 0; + } + if(block->header_contents) + { + free(block->header_contents); + block->header_contents = 0; + } + if(block->block_contents) + { + free(block->block_contents); + block->block_contents = 0; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_atom_name(struct trg_trajectory *trg_data, + struct trg_atom *atom, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(atom->name && strlen(atom->name) < len) + { + free(atom->name); + atom->name = 0; + } + if(!atom->name) + { + atom->name = (char *) malloc(len); + if(!atom->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(atom->name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_atom_type(struct trg_trajectory *trg_data, + struct trg_atom *atom, + const char *new_type) +{ + int len; + + len = min(strlen(new_type) + 1, TRG_MAX_STR_LEN); + + if(atom->atom_type && strlen(atom->atom_type) < len) + { + free(atom->atom_type); + atom->atom_type = 0; + } + if(!atom->atom_type) + { + atom->atom_type = (char *) malloc(len); + if(!atom->atom_type) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(atom->atom_type, new_type, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_init_atom(struct trg_atom *atom) +{ + atom->name = 0; + atom->atom_type = 0; + + return(TRG_SUCCESS); +} + +trg_function_status trg_destroy_atom(struct trg_atom *atom) +{ + if(atom->name) + { + free(atom->name); + atom->name = 0; + } + if(atom->atom_type) + { + free(atom->atom_type); + atom->atom_type = 0; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_add_molecule(struct trg_trajectory *trg_data, + const char *name, + struct trg_molecule **molecule) +{ + struct trg_molecule *new_molecules; + uint64_t *new_molecule_cnt_list; + int id, i; + trg_bool found_id = TRUE; + + new_molecules = (struct trg_molecule *)realloc(trg_data->molecules, + sizeof(struct trg_molecule) * + (trg_data->n_molecules + 1)); + + if(!new_molecules) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(struct trg_molecule) * (trg_data->n_molecules + 1), + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + new_molecule_cnt_list = (uint64_t *) realloc(trg_data->molecule_cnt_list, + sizeof(uint64_t) * + (trg_data->n_molecules + 1)); + + if(!new_molecules) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(uint64_t) * (trg_data->n_molecules + 1), + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + trg_data->molecules = new_molecules; + trg_data->molecule_cnt_list = new_molecule_cnt_list; + + *molecule = &new_molecules[trg_data->n_molecules]; + + trg_init_molecule(*molecule); + trg_set_molecule_name(trg_data, *molecule, name); + + /* FIXME: Should this be a function argument instead? */ + trg_data->molecule_cnt_list[trg_data->n_molecules] = 0; + + /* Find an unused ID */ + id = 0; + while(found_id) + { + found_id = FALSE; + for(i = trg_data->n_molecules; i--;) + { + if(trg_data->molecules[i].id == id) + { + found_id = TRUE; + i = 0; + } + } + if(found_id) + { + id++; + } + } + + (*molecule)->id = id; + + trg_data->n_molecules++; + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_molecule_name(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(molecule->name && strlen(molecule->name) < len) + { + free(molecule->name); + molecule->name = 0; + } + if(!molecule->name) + { + molecule->name = (char *) malloc(len); + if(!molecule->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(molecule->name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_get_molecule_cnt(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + uint64_t *cnt) +{ + int i, index = -1; + + for(i = trg_data->n_molecules; i--;) + { + if(&trg_data->molecules[i] == molecule) + { + index = i; + i = 0; + } + } + if(index == -1) + { + return(TRG_FAILURE); + } + *cnt = trg_data->molecule_cnt_list[index]; + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_molecule_cnt(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + const uint64_t cnt) +{ + int i, index = -1, old_cnt; + + for(i = trg_data->n_molecules; i--;) + { + if(&trg_data->molecules[i] == molecule) + { + index = i; + i = 0; + } + } + if(index == -1) + { + return(TRG_FAILURE); + } + old_cnt = trg_data->molecule_cnt_list[index]; + trg_data->molecule_cnt_list[index] = cnt; + + trg_data->n_particles += (cnt-old_cnt) * + trg_data->molecules[index].n_atoms; + + return(TRG_SUCCESS); +} + +trg_function_status trg_add_chain_to_molecule(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + const char *name, + struct trg_chain **chain) +{ + struct trg_chain *new_chains; + + new_chains = (struct trg_chain *) realloc(molecule->chains, + sizeof(struct trg_chain) * + (molecule->n_chains + 1)); + + if(!new_chains) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(struct trg_chain) * (molecule->n_chains + 1), + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + molecule->chains = new_chains; + + *chain = &new_chains[molecule->n_chains]; + (*chain)->name = 0; + + trg_set_chain_name(trg_data, *chain, name); + + (*chain)->molecule = molecule; + (*chain)->id = molecule->n_chains; + (*chain)->n_residues = 0; + + molecule->n_chains++; + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_chain_name(struct trg_trajectory *trg_data, + struct trg_chain *chain, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(chain->name && strlen(chain->name) < len) + { + free(chain->name); + chain->name = 0; + } + if(!chain->name) + { + chain->name = (char *) malloc(len); + if(!chain->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(chain->name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_add_residue_to_chain(struct trg_trajectory *trg_data, + struct trg_chain *chain, + const char *name, + struct trg_residue **residue) +{ + int curr_index; + struct trg_residue *new_residues, *temp_residue, *last_residue; + struct trg_molecule *molecule = chain->molecule; + + if(chain->n_residues) + { + curr_index = (chain->residues - molecule->residues) / + sizeof(struct trg_residue); + } + else + { + curr_index = -1; + } + + new_residues = (struct trg_residue *) realloc(molecule->residues, + sizeof(struct trg_residue) * + (molecule->n_residues + 1)); + + if(!new_residues) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(struct trg_residue) * (molecule->n_residues + 1), + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + molecule->residues = new_residues; + + if(curr_index != -1) + { + chain->residues = new_residues + curr_index * sizeof(struct trg_residue); + if(molecule->n_residues) + { + last_residue = &new_residues[molecule->n_atoms-1]; + + temp_residue = chain->residues + (chain->n_residues - 1); + /* Make space in list of residues to add the new residues together with the other + * residues of this chain */ + if(temp_residue != last_residue) + { + temp_residue++; + memmove(temp_residue + 1, temp_residue, + last_residue - temp_residue); + } + } + } + else + { + curr_index = molecule->n_residues; + } + + *residue = &new_residues[curr_index + chain->n_residues]; + + if(!chain->n_residues) + { + chain->residues = *residue; + } + + (*residue)->name = 0; + trg_set_residue_name(trg_data, *residue, name); + + (*residue)->chain = chain; + (*residue)->id = chain->n_residues; + (*residue)->n_atoms = 0; + (*residue)->atoms = 0; + + chain->n_residues++; + molecule->n_residues++; + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_residue_name(struct trg_trajectory *trg_data, + struct trg_residue *residue, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(residue->name && strlen(residue->name) < len) + { + free(residue->name); + residue->name = 0; + } + if(!residue->name) + { + residue->name = (char *) malloc(len); + if(!residue->name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(residue->name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_add_atom_to_residue(struct trg_trajectory *trg_data, + struct trg_residue *residue, + const char *atom_name, + const char *atom_type, + struct trg_atom **atom) +{ + int curr_index; + struct trg_atom *new_atoms, *temp_atom, *last_atom; + struct trg_molecule *molecule = residue->chain->molecule; + + if(residue->n_atoms) + { + curr_index = (residue->atoms - molecule->atoms) / + sizeof(struct trg_atom); + } + else + { + curr_index = -1; + } + + new_atoms = (struct trg_atom *) realloc(molecule->atoms, + sizeof(struct trg_atom) * + (molecule->n_atoms + 1)); + + if(!new_atoms) + { + printf("Cannot allocate memory (%"PRId64" bytes). %s: %d\n", + sizeof(struct trg_atom) * (molecule->n_atoms + 1), + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + molecule->atoms = new_atoms; + + if(curr_index != -1) + { + residue->atoms = new_atoms + curr_index * sizeof(struct trg_atom); + if(molecule->n_atoms) + { + last_atom = &new_atoms[molecule->n_atoms-1]; + + temp_atom = residue->atoms + (residue->n_atoms - 1); + /* Make space in list of atoms to add the new atoms together with the other + * atoms of this residue */ + if(temp_atom != last_atom) + { + temp_atom++; + memmove(temp_atom + 1, temp_atom, + last_atom - temp_atom); + } + } + } + else + { + curr_index = molecule->n_atoms; + } + + *atom = &new_atoms[curr_index + residue->n_atoms]; + + if(!residue->n_atoms) + { + residue->atoms = *atom; + } + + trg_init_atom(*atom); + trg_set_atom_name(trg_data, *atom, atom_name); + trg_set_atom_type(trg_data, *atom, atom_type); + + (*atom)->residue = residue; + (*atom)->id = molecule->n_atoms; + + residue->n_atoms++; + molecule->n_atoms++; + + return(TRG_SUCCESS); +} + + +trg_function_status trg_init_molecule(struct trg_molecule *molecule) +{ + molecule->name = 0; + molecule->n_chains = 0; + molecule->chains = 0; + molecule->n_residues = 0; + molecule->residues = 0; + molecule->n_atoms = 0; + molecule->atoms = 0; + molecule->n_bonds = 0; + molecule->bonds = 0; + + return(TRG_SUCCESS); +} + +trg_function_status trg_destroy_molecule(struct trg_molecule *molecule) +{ + int i; + + if(molecule->name) + { + free(molecule->name); + molecule->name = 0; + } + + if(molecule->chains) + { + for(i = molecule->n_chains; i--;) + { + if(molecule->chains[i].name) + { + free(molecule->chains[i].name); + molecule->chains[i].name = 0; + } + } + free(molecule->chains); + molecule->chains = 0; + } + molecule->n_chains = 0; + + if(molecule->residues) + { + for(i = molecule->n_residues; i--;) + { + if(molecule->residues[i].name) + { + free(molecule->residues[i].name); + molecule->residues[i].name = 0; + } + } + free(molecule->residues); + molecule->residues = 0; + } + molecule->n_residues = 0; + + if(molecule->atoms) + { + for(i = molecule->n_atoms; i--;) + { + trg_destroy_atom(&molecule->atoms[i]); + } + free(molecule->atoms); + molecule->atoms = 0; + } + molecule->n_atoms = 0; + + if(molecule->bonds) + { + free(molecule->bonds); + molecule->bonds = 0; + } + molecule->n_bonds = 0; + + return(TRG_SUCCESS); +} + +trg_function_status trg_init_trajectory(struct trg_trajectory *trg_data) +{ + time_t seconds; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + trg_data->input_file_path = 0; + trg_data->input_file = 0; + trg_data->input_file_pos = 0; + trg_data->input_file_len = 0; + trg_data->output_file_path = 0; + trg_data->output_file = 0; + trg_data->output_file_pos = 0; + + trg_data->program_name = 0; + trg_data->forcefield_name = 0; + + /* FIXME: No unistd.h on Windows!! */ +// trg_data->user_name = (char *) malloc(LOGIN_NAME_MAX); +// if(getlogin_r(trg_data->user_name, LOGIN_NAME_MAX) != 0) +// { +// printf("Cannot get user name. %s: %d\n", __FILE__, __LINE__); +// free(trg_data->user_name); +// trg_data->user_name = 0; +// } + trg_data->user_name = 0; + + seconds = time(0); + if ( seconds == -1) + { + printf("Cannot get time. %s: %d\n", __FILE__, __LINE__); + } + else + { + trg_data->time = seconds; + } + + /* FIXME: No unistd.h on Windows!! */ + /* FIXME: Append operating system to computer_name */ +// trg_data->computer_name = (char *) malloc(HOST_NAME_MAX); +// if(gethostname(trg_data->computer_name, HOST_NAME_MAX) != 0) +// { +// printf("Cannot get computer name. %s: %d\n", __FILE__, __LINE__); +// free(trg_data->computer_name); +// trg_data->computer_name = 0; +// } + trg_data->computer_name = 0; + + trg_data->pgp_signature = 0; + trg_data->var_num_atoms_flag = TRG_CONSTANT_N_ATOMS; + trg_data->first_trajectory_frame_set_input_file_pos = -1; + trg_data->last_trajectory_frame_set_input_file_pos = -1; + trg_data->current_trajectory_frame_set_input_file_pos = -1; + trg_data->first_trajectory_frame_set_output_file_pos = -1; + trg_data->last_trajectory_frame_set_output_file_pos = -1; + trg_data->current_trajectory_frame_set_output_file_pos = -1; + trg_data->frame_set_n_frames = 100; + trg_data->n_trajectory_frame_sets = 0; + trg_data->n_trajectory_blocks = 0; + trg_data->stride_length = 100; + + trg_data->n_particle_data_blocks = 0; + trg_data->n_data_blocks = 0; + + trg_data->non_tr_particle_data = 0; + trg_data->non_tr_data = 0; + + frame_set->contents.n_blocks = 0; + frame_set->contents.block_names = 0; + frame_set->n_mapping_blocks = 0; + frame_set->mappings = 0; + frame_set->molecule_cnt_list = 0; + + frame_set->n_particle_data_blocks = 0; + frame_set->n_data_blocks = 0; + + frame_set->tr_particle_data = 0; + frame_set->tr_data = 0; + + trg_data->n_molecules = 0; + trg_data->molecules = 0; + trg_data->molecule_cnt_list = 0; + trg_data->n_particles = 0; + + trg_data->n_id_name_pairs = 0; + trg_data->id_name_pairs = 0; + + /* Check the endianness of the computer */ + static uint32_t endianness_32 = 0x01234567; + /* 0x01234567 */ + if ( *(const uint8_t*)&endianness_32 == 0x01 ) + { + trg_data->endianness_32 = TRG_BIG_ENDIAN_32; + } + + /* 0x67452301 */ + else if( *(const uint8_t*)&endianness_32 == 0x67 ) + { + trg_data->endianness_32 = TRG_LITTLE_ENDIAN_32; + + } + + /* 0x45670123 */ + else if ( *(const uint8_t*)&endianness_32 == 0x45 ) + { + trg_data->endianness_32 = TRG_BYTE_PAIR_SWAP_32; + } + + static uint64_t endianness_64 = 0x0123456789ABCDEF; + /* 0x0123456789ABCDEF */ + if ( *(const uint8_t*)&endianness_64 == 0x01 ) + { + trg_data->endianness_64 = TRG_BIG_ENDIAN_64; + } + + /* 0xEFCDAB8967452301 */ + else if ( *(const uint8_t*)&endianness_64 == 0xEF ) + { + trg_data->endianness_64 = TRG_LITTLE_ENDIAN_64; + } + + /* 0x89ABCDEF01234567 */ + else if ( *(const uint8_t*)&endianness_64 == 0x89 ) + { + trg_data->endianness_64 = TRG_QUAD_SWAP_64; + } + + /* 0x45670123CDEF89AB */ + else if ( *(const uint8_t*)&endianness_64 == 0x45 ) + { + trg_data->endianness_64 = TRG_BYTE_PAIR_SWAP_64; + } + + /* 0x23016745AB89EFCD */ + else if ( *(const uint8_t*)&endianness_64 == 0x23 ) + { + trg_data->endianness_64 = TRG_BYTE_SWAP_64; + } + + + + trg_init_block(&trg_data->non_trajectory_blocks[0]); + trg_data->non_trajectory_blocks[0].id = TRG_GENERAL_INFO; + trg_set_block_name(trg_data, &trg_data->non_trajectory_blocks[0], + "GENERAL INFO"); + + trg_data->current_trajectory_frame_set.next_frame_set_file_pos = -1; + trg_data->current_trajectory_frame_set.prev_frame_set_file_pos = -1; + + /* The Endianness and String Length block and the Trajectory Info block + * are present. */ + trg_data->n_non_trajectory_blocks = 1; + + return(TRG_SUCCESS); +} + +trg_function_status trg_destroy_trajectory(struct trg_trajectory *trg_data) +{ + uint64_t n_frames, n_particles; + int i, j, k, l; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + struct trg_particle_mapping *mapping; + + if(trg_data->input_file_path) + { + free(trg_data->input_file_path); + trg_data->input_file_path = 0; + } + + if(trg_data->input_file) + { + fclose(trg_data->input_file); + trg_data->input_file = 0; + } + + if(trg_data->output_file_path) + { + free(trg_data->output_file_path); + trg_data->output_file_path = 0; + } + + if(trg_data->output_file) + { + fclose(trg_data->output_file); + trg_data->output_file = 0; + } + + if(trg_data->program_name) + { + free(trg_data->program_name); + trg_data->program_name = 0; + } + + if(trg_data->forcefield_name) + { + free(trg_data->forcefield_name); + trg_data->forcefield_name = 0; + } + + if(trg_data->user_name) + { + free(trg_data->user_name); + trg_data->user_name = 0; + } + + if(trg_data->computer_name) + { + free(trg_data->computer_name); + trg_data->computer_name = 0; + } + + if(trg_data->pgp_signature) + { + free(trg_data->pgp_signature); + trg_data->pgp_signature = 0; + } + + if(frame_set->contents.block_names) + { + for(i = frame_set->contents.n_blocks; i--;) + { + if(frame_set->contents.block_names[i]) + { + free(frame_set->contents.block_names[i]); + frame_set->contents.block_names[i] = 0; + } + } + free(frame_set->contents.block_names); + frame_set->contents.block_names = 0; + } + + for(i = frame_set->n_mapping_blocks; i--;) + { + mapping = &frame_set->mappings[i]; + if(mapping->real_particle_numbers) + { + free(mapping->real_particle_numbers); + mapping->real_particle_numbers = 0; + } + } + free(frame_set->mappings); + frame_set->n_mapping_blocks = 0; + + if(frame_set->molecule_cnt_list) + { + free(frame_set->molecule_cnt_list); + frame_set->molecule_cnt_list = 0; + } + + for(i=trg_data->n_non_trajectory_blocks; i--;) + { + trg_destroy_block(&trg_data->non_trajectory_blocks[i]); + } + trg_data->n_trajectory_blocks = 0; + + if(trg_data->var_num_atoms_flag) + { + n_particles = trg_data->current_trajectory_frame_set.n_particles; + } + else + { + n_particles = trg_data->n_particles; + } + + if(trg_data->non_tr_particle_data) + { + for(i = trg_data->n_particle_data_blocks; i--; ) + { + if(trg_data->non_tr_particle_data[i].values) + { + /* Only one frame for non-trajectory data */ + j = 0; + if(trg_data->non_tr_particle_data[i].values[j]) + { + for(k = n_particles; k--;) + { + if(trg_data->non_tr_particle_data[i].values[j][k]) + { + if(trg_data->non_tr_particle_data[i].datatype == + TRG_CHAR_DATA) + { + for(l = trg_data->non_tr_particle_data[i]. + n_values_per_frame; + l--;) + { + if(trg_data->non_tr_particle_data[i]. + values[j][k][l].c) + { + free(trg_data->non_tr_particle_data[i]. + values[j][k][l].c); + trg_data->non_tr_particle_data[i]. + values[j][k][l].c = 0; + } + } + } + free(trg_data->non_tr_particle_data[i]. + values[j][k]); + trg_data->non_tr_particle_data[i]. + values[j][k] = 0; + } + } + free(trg_data->non_tr_particle_data[i].values[j]); + trg_data->non_tr_particle_data[i].values[j] = 0; + } + free(trg_data->non_tr_particle_data[i].values); + trg_data->non_tr_particle_data[i].values = 0; + } + if(trg_data->non_tr_particle_data[i].block_name) + { + free(trg_data->non_tr_particle_data[i].block_name); + trg_data->non_tr_particle_data[i].block_name = 0; + } + } + free(trg_data->non_tr_particle_data); + trg_data->non_tr_particle_data = 0; + } + + if(trg_data->non_tr_data) + { + for(i = trg_data->n_data_blocks; i--;) + { + if(trg_data->non_tr_data[i].values) + { + /* Only one frame for non-trajectory data */ + if(trg_data->non_tr_data[i].values[0]) + { + if(trg_data->non_tr_data[i].datatype == + TRG_CHAR_DATA) + { + for(k = trg_data->non_tr_data[i].n_values_per_frame; + k--;) + { + if(trg_data->non_tr_data[i].values[0][k].c) + { + free(trg_data->non_tr_data[i].values[0][k].c); + trg_data->non_tr_data[i].values[0][k].c = 0; + } + } + } + free(trg_data->non_tr_data[i].values[0]); + trg_data->non_tr_data[i].values[0] = 0; + } + free(trg_data->non_tr_data[i].values); + trg_data->non_tr_data[i].values = 0; + } + if(trg_data->non_tr_data[i].block_name) + { + free(trg_data->non_tr_data[i].block_name); + trg_data->non_tr_data[i].block_name = 0; + } + } + free(trg_data->non_tr_data); + trg_data->non_tr_data = 0; + } + + trg_data->n_particle_data_blocks = 0; + trg_data->n_data_blocks = 0; + + if(frame_set->tr_particle_data) + { + for(i = frame_set->n_particle_data_blocks; i--; ) + { + if(frame_set->tr_particle_data[i].values) + { + n_frames = max(1, frame_set->tr_particle_data[i].n_frames); + for(j = n_frames; j--;) + { + if(frame_set->tr_particle_data[i].values[j]) + { + for(k = n_particles; k--;) + { + if(frame_set->tr_particle_data[i].values[j][k]) + { + if(frame_set->tr_particle_data[i].datatype == + TRG_CHAR_DATA) + { + for(l = frame_set->tr_particle_data[i]. + n_values_per_frame; + l--;) + { + if(frame_set->tr_particle_data[i]. + values[j][k][l].c) + { + free(frame_set->tr_particle_data[i]. + values[j][k][l].c); + + frame_set->tr_particle_data[i]. + values[j][k][l].c = 0; + } + } + } + free(frame_set->tr_particle_data[i]. + values[j][k]); + frame_set->tr_particle_data[i].values[j][k] = 0; + } + } + free(frame_set->tr_particle_data[i].values[j]); + frame_set->tr_particle_data[i].values[j] = 0; + } + } + free(frame_set->tr_particle_data[i].values); + frame_set->tr_particle_data[i].values = 0; + } + if(frame_set->tr_particle_data[i].block_name) + { + free(frame_set->tr_particle_data[i].block_name); + frame_set->tr_particle_data[i].block_name = 0; + } + } + free(frame_set->tr_particle_data); + frame_set->tr_particle_data = 0; + } + + if(frame_set->tr_data) + { + for(i = frame_set->n_data_blocks; i--;) + { + if(frame_set->tr_data[i].values) + { + n_frames = max(1, frame_set->tr_data[i].n_frames); + for(j = n_frames; j--;) + { + if(frame_set->tr_data[i].values[j]) + { + if(frame_set->tr_data[i].datatype == + TRG_CHAR_DATA) + { + for(k = frame_set->tr_data[i].n_values_per_frame; + k--;) + { + if(frame_set->tr_data[i].values[j][k].c) + { + free(frame_set->tr_data[i]. + values[j][k].c); + + frame_set->tr_data[i]. + values[j][k].c = 0; + } + } + } + free(frame_set->tr_data[i].values[j]); + frame_set->tr_data[i].values[j] = 0; + } + } + free(frame_set->tr_data[i].values); + frame_set->tr_data[i].values = 0; + } + if(frame_set->tr_data[i].block_name) + { + free(frame_set->tr_data[i].block_name); + frame_set->tr_data[i].block_name = 0; + } + } + free(frame_set->tr_data); + frame_set->tr_data = 0; + } + + frame_set->n_particle_data_blocks = 0; + frame_set->n_data_blocks = 0; + + if(trg_data->molecules) + { + for(i=trg_data->n_molecules; i--;) + { + trg_destroy_molecule(&trg_data->molecules[i]); + } + free(trg_data->molecules); + trg_data->molecules = 0; + trg_data->n_molecules = 0; + } + if(trg_data->molecule_cnt_list) + { + free(trg_data->molecule_cnt_list); + trg_data->molecule_cnt_list = 0; + } + return(TRG_SUCCESS); +} + +trg_function_status trg_set_input_file(struct trg_trajectory *trg_data, + const char *file_name) +{ + int len; + char *temp; + + if(trg_data->input_file_path && strcmp(trg_data->input_file_path, + file_name) == 0) + { + return(TRG_SUCCESS); + } + + if(trg_data->input_file) + { + fclose(trg_data->input_file); + } + + len = min(strlen(file_name) + 1, TRG_MAX_STR_LEN); + temp = realloc(trg_data->input_file_path, len); + if(!temp) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + trg_data->input_file_path = temp; + + strncpy(trg_data->input_file_path, file_name, len); + + return(trg_init_input_file(trg_data, FALSE)); +} + +trg_function_status trg_set_output_file(struct trg_trajectory *trg_data, + const char *file_name) +{ + int len; + char *temp; + + if(trg_data->output_file_path && + strcmp(trg_data->output_file_path, file_name) == 0) + { + return(TRG_SUCCESS); + } + + if(trg_data->output_file) + { + fclose(trg_data->output_file); + } + + len = min(strlen(file_name) + 1, TRG_MAX_STR_LEN); + temp = realloc(trg_data->output_file_path, len); + if(!temp) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + trg_data->output_file_path = temp; + + strncpy(trg_data->output_file_path, file_name, len); + + return(trg_init_output_file(trg_data, FALSE)); +} + + +trg_function_status trg_set_program_name(struct trg_trajectory *trg_data, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(trg_data->program_name && strlen(trg_data->program_name) < len) + { + free(trg_data->program_name); + trg_data->program_name = 0; + } + if(!trg_data->program_name) + { + trg_data->program_name = (char *) malloc(len); + if(!trg_data->program_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(trg_data->program_name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_forcefield_name(struct trg_trajectory *trg_data, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(trg_data->forcefield_name && strlen(trg_data->forcefield_name) < len) + { + free(trg_data->forcefield_name); + trg_data->forcefield_name = 0; + } + if(!trg_data->forcefield_name) + { + trg_data->forcefield_name = (char *) malloc(len); + if(!trg_data->forcefield_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(trg_data->forcefield_name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_user_name(struct trg_trajectory *trg_data, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(trg_data->user_name && strlen(trg_data->user_name) < len) + { + free(trg_data->user_name); + trg_data->user_name = 0; + } + if(!trg_data->user_name) + { + trg_data->user_name = (char *) malloc(len); + if(!trg_data->user_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(trg_data->user_name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_computer_name(struct trg_trajectory *trg_data, + const char *new_name) +{ + int len; + + len = min(strlen(new_name) + 1, TRG_MAX_STR_LEN); + + if(trg_data->computer_name && strlen(trg_data->computer_name) < len) + { + free(trg_data->computer_name); + trg_data->computer_name = 0; + } + if(!trg_data->computer_name) + { + trg_data->computer_name = (char *) malloc(len); + if(!trg_data->computer_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(trg_data->computer_name, new_name, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_set_signature(struct trg_trajectory *trg_data, + const char *signature) +{ + int len; + + len = min(strlen(signature) + 1, TRG_MAX_STR_LEN); + + if(trg_data->pgp_signature && strlen(trg_data->pgp_signature) < len) + { + free(trg_data->pgp_signature); + trg_data->pgp_signature = 0; + } + if(!trg_data->pgp_signature) + { + trg_data->pgp_signature = (char *) malloc(len); + if(!trg_data->pgp_signature) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + strncpy(trg_data->pgp_signature, signature, len); + + return(TRG_SUCCESS); +} + +trg_function_status trg_read_file_headers(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + int i, cnt = 0, prev_pos = 0; + struct trg_gen_block *block = trg_data->non_trajectory_blocks; + + trg_data->input_file_pos = 0; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + if(!trg_data->input_file_len) + { + fseek(trg_data->input_file, 0, SEEK_END); + trg_data->input_file_len = ftell(trg_data->input_file); + fseek(trg_data->input_file, 0, SEEK_SET); + } + + for(i = trg_data->n_non_trajectory_blocks; i--;) + { + trg_destroy_block(block++); + } + trg_data->n_non_trajectory_blocks = 0; + + block = trg_data->non_trajectory_blocks; + + trg_init_block(block); + /* Non trajectory blocks (they come before the trajectory + * blocks in the file) */ + while (prev_pos < trg_data->input_file_len && + trg_read_block_header(trg_data, block) != TRG_CRITICAL && + block->id != -1 && + block->id != TRG_TRAJECTORY_FRAME_SET && + trg_data->n_non_trajectory_blocks < 32) + { +// printf("Reading block header %d: %s\n", (int)block->id, block->name); + if(trg_read_next_block(trg_data, block, + TRG_KEEP_FILE_OPEN) == TRG_SUCCESS) + { +// printf("Read block %s\n", block->name); + block++; + cnt++; + trg_data->n_non_trajectory_blocks++; + if(trg_init_block(block) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + } + else + { + trg_destroy_block(block); + } + prev_pos = ftell(trg_data->input_file); + } + + /* Go back if a trajectory block was encountered */ + if(block->id == TRG_TRAJECTORY_FRAME_SET) + { + trg_destroy_block(block); + fseek(trg_data->input_file, prev_pos, SEEK_SET); + } + + if(close_file) + { + trg_data->input_file_pos=ftell(trg_data->input_file); + fclose(trg_data->input_file); + trg_data->input_file = 0; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_write_file_headers(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + int i; + struct trg_gen_block *block, data_block; + + trg_data->output_file_pos = 0; + + if(trg_init_output_file(trg_data, TRUE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + + for(i=0; i<trg_data->n_non_trajectory_blocks; ++i) + { + block = &trg_data->non_trajectory_blocks[i]; + if(block->id == TRG_GENERAL_INFO) + { + if(trg_write_general_info_block(trg_data, block, TRG_NORMAL_WRITE) + != TRG_SUCCESS) + { + printf("Error writing general info block of file %s. %s: %d\n", + trg_data->input_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + break; + } + } + + for(i=0; i<trg_data->n_non_trajectory_blocks; ++i) + { + block = &trg_data->non_trajectory_blocks[i]; + if(block->id == TRG_MOLECULES) + { + if(trg_write_molecules_block(trg_data, block, TRG_NORMAL_WRITE) + != TRG_SUCCESS) + { + printf("Error writing atom names block of file %s. %s: %d\n", + trg_data->input_file_path, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + break; + } + } + + /* FIXME: Currently writing non-trajectory data blocks here. + * Should perhaps be moved. */ + trg_init_block(&data_block); + for(i = 0; i < trg_data->n_data_blocks; i++) + { + data_block.id = trg_data->non_tr_data[i].block_id; + trg_write_data_block(trg_data, &data_block, + i, TRG_NORMAL_WRITE); + } + + for(i = 0; i < trg_data->n_particle_data_blocks; i++) + { + data_block.id = trg_data->non_tr_particle_data[i].block_id; + trg_write_particle_data_block(trg_data, &data_block, + i, 0, TRG_NORMAL_WRITE); + } + + trg_destroy_block(&data_block); + + if(close_file) + { + trg_data->output_file_pos=ftell(trg_data->output_file); + fclose(trg_data->output_file); + trg_data->output_file = 0; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_read_next_block(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + trg_close_file_flag close_file) +{ + switch(block->id) + { + case TRG_TRAJECTORY_FRAME_SET: + return(trg_read_frame_set_block(trg_data, block)); + case TRG_BLOCK_TABLE_OF_CONTENTS: + return(trg_read_trajectory_toc_block(trg_data, block)); + case TRG_PARTICLE_MAPPING: + return(trg_read_trajectory_mapping_block(trg_data, block)); + case TRG_GENERAL_INFO: + return(trg_read_general_info_block(trg_data, block)); + case TRG_MOLECULES: + return(trg_read_molecules_block(trg_data, block)); + default: + if(block->id >= TRG_TRAJ_BOX_SHAPE) + { + return(trg_read_data_block_contents(trg_data, block)); + } + else + { + /* Skip to the next block */ + fseek(trg_data->input_file, block->block_contents_size, SEEK_CUR); + return(TRG_FAILURE); + } + } + + /* FIXME: Never reached. */ + if(close_file) + { + trg_data->input_file_pos=ftell(trg_data->input_file); + fclose(trg_data->input_file); + trg_data->input_file = 0; + } +} + +// trg_function_status trg_write_block(struct trg_trajectory *trg_data, +// struct trg_gen_block *block, +// trg_close_file_flag close_file) +// { +// if(trg_data->output_file) +// { +// trg_data->output_file_pos = ftell(trg_data->output_file); +// } +// +// switch(block->id) +// { +// case TRG_TRAJECTORY_FRAME_SET: +// return(trg_write_frame_set_block(trg_data, block, TRG_NORMAL_WRITE)); +// break; +// case TRG_BLOCK_TABLE_OF_CONTENTS: +// return(trg_write_trajectory_toc_block(trg_data, block, +// TRG_NORMAL_WRITE)); +// break; +// case TRG_PARTICLE_MAPPING: +// return(trg_write_trajectory_mapping_block(trg_data, block, +// TRG_NORMAL_WRITE)); +// break; +// case TRG_TRAJ_POSITIONS: +// break; +// case TRG_TRAJ_VELOCITIES: +// break; +// case TRG_TRAJ_FORCES: +// break; +// case TRG_TRAJ_BOX_SHAPE: +// break; +// case TRG_GENERAL_INFO: +// return(trg_write_general_info_block(trg_data, block, +// TRG_NORMAL_WRITE)); +// case TRG_MOLECULES: +// return(trg_write_molecules_block(trg_data, block, +// TRG_NORMAL_WRITE)); +// default: +// if(block->id > TRG_TRAJ_FORCES) +// { +// /* Implement writing data blocks. */ +// return(TRG_FAILURE); +// } +// else +// { +// return(TRG_FAILURE); +// } +// } +// +// /* FIXME: Never reached. */ +// if(close_file) +// { +// fclose(trg_data->output_file); +// trg_data->output_file = 0; +// } +// +// return(TRG_SUCCESS); +// } + +trg_function_status trg_read_next_frame_set(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + long int file_pos; + struct trg_gen_block block; + trg_function_status stat = TRG_SUCCESS; + + if(trg_init_input_file(trg_data, FALSE) != TRG_SUCCESS) + { + return(TRG_CRITICAL); + } + + trg_init_block(&block); + + file_pos = trg_data->current_trajectory_frame_set.next_frame_set_file_pos; + + if(file_pos > 0) + { + fseek(trg_data->input_file, + file_pos, + SEEK_SET); + } + else + { + return(TRG_CRITICAL); + } + + if(!trg_data->input_file_len) + { + fseek(trg_data->input_file, 0, SEEK_END); + trg_data->input_file_len = ftell(trg_data->input_file); + fseek(trg_data->input_file, file_pos, SEEK_SET); + } + + stat = trg_read_block_header(trg_data, &block); + if(stat == TRG_CRITICAL || block.id != TRG_TRAJECTORY_FRAME_SET) + { + return(TRG_CRITICAL); + } + + trg_data->current_trajectory_frame_set_input_file_pos = file_pos; + + if(trg_read_next_block(trg_data, &block, + TRG_KEEP_FILE_OPEN) == TRG_SUCCESS) + { + file_pos = ftell(trg_data->input_file); + /* Read all blocks until next frame set block */ + stat = trg_read_block_header(trg_data, &block); + while(file_pos < trg_data->input_file_len && + stat != TRG_CRITICAL && + block.id != TRG_TRAJECTORY_FRAME_SET) + { + stat = trg_read_next_block(trg_data, &block, + TRG_KEEP_FILE_OPEN) == TRG_SUCCESS; + + if(stat != TRG_CRITICAL) + { + file_pos = ftell(trg_data->input_file); + if(file_pos < trg_data->input_file_len) + { + stat = trg_read_block_header(trg_data, &block); + } + } + } + if(stat == TRG_CRITICAL) + { + return(stat); + } + + if(block.id == TRG_TRAJECTORY_FRAME_SET) + { + fseek(trg_data->input_file, file_pos, SEEK_SET); + } + } + + trg_data->input_file_pos=ftell(trg_data->input_file); + + trg_destroy_block(&block); + + if(close_file) + { + fclose(trg_data->input_file); + trg_data->input_file = 0; + } + + return(TRG_SUCCESS); +} + +trg_function_status trg_write_frame_set(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + int i, j; + struct trg_gen_block block; + struct trg_trajectory_frame_set *frame_set = + &trg_data->current_trajectory_frame_set; + + trg_function_status stat; + + + if(trg_data->output_file) + { + trg_data->current_trajectory_frame_set_output_file_pos = + ftell(trg_data->output_file); + } + else + { + trg_data->current_trajectory_frame_set_output_file_pos = + trg_data->output_file_pos; + } + + trg_init_block(&block); + block.id = TRG_TRAJECTORY_FRAME_SET; + + trg_write_frame_set_block(trg_data, &block, TRG_NORMAL_WRITE); + + if(frame_set->contents.n_blocks > 0) + { + block.id = TRG_BLOCK_TABLE_OF_CONTENTS; + trg_write_trajectory_toc_block(trg_data, &block, TRG_NORMAL_WRITE); + } + for(i = 0; i<frame_set->n_data_blocks; i++) + { + block.id = frame_set->tr_data[i].block_id; + trg_write_data_block(trg_data, &block, i, TRG_NORMAL_WRITE); + } + if(frame_set->n_mapping_blocks) + { + for(i = 0; i < frame_set->n_mapping_blocks; i++) + { + block.id = TRG_PARTICLE_MAPPING; + if(frame_set->mappings[i].n_particles > 0) + { + trg_write_trajectory_mapping_block(trg_data, &block, i, + TRG_NORMAL_WRITE); + for(j = 0; j<frame_set->n_particle_data_blocks; j++) + { + block.id = frame_set->tr_particle_data[i].block_id; + trg_write_particle_data_block(trg_data, &block, + j, &frame_set->mappings[i], + TRG_NORMAL_WRITE); + } + } + } + } + else + { + for(i = 0; i<frame_set->n_particle_data_blocks; i++) + { + block.id = frame_set->tr_particle_data[i].block_id; + trg_write_particle_data_block(trg_data, &block, + i, 0, TRG_NORMAL_WRITE); + } + } + + trg_data->output_file_pos = ftell(trg_data->output_file); + + stat = trg_update_header_pointers(trg_data); + + if(stat == TRG_SUCCESS) + { + stat = trg_update_frame_set_pointers(trg_data); + } + + trg_destroy_block(&block); + + if(close_file) + { + fclose(trg_data->input_file); + trg_data->input_file = 0; + } + + return(stat); +} + +trg_function_status trg_new_frame_set(struct trg_trajectory *trg_data, + const uint64_t first_frame, + const uint64_t n_frames) +{ + int i; + struct trg_gen_block block; + struct trg_trajectory_frame_set *frame_set; + struct trg_particle_mapping *mapping; + FILE *temp = trg_data->input_file; + uint64_t curr_pos; + + frame_set = &trg_data->current_trajectory_frame_set; + + if(trg_data->n_trajectory_frame_sets) + { + frame_set->prev_frame_set_file_pos = + trg_data->current_trajectory_frame_set_output_file_pos; + } + + if(trg_data->output_file) + { + trg_data->current_trajectory_frame_set_output_file_pos = + ftell(trg_data->output_file); + } + else + { + trg_data->current_trajectory_frame_set_output_file_pos = + trg_data->output_file_pos; + } + + if(frame_set->n_mapping_blocks && frame_set->mappings) + { + for(i = frame_set->n_mapping_blocks; i--;) + { + mapping = &frame_set->mappings[i]; + if(mapping->real_particle_numbers) + { + free(mapping->real_particle_numbers); + mapping->real_particle_numbers = 0; + } + } + free(frame_set->mappings); + frame_set->n_mapping_blocks = 0; + } + + trg_data->n_trajectory_frame_sets++; + + /* Set the long range pointers */ + if(trg_data->n_trajectory_frame_sets == trg_data->stride_length + 1) + { + frame_set->long_stride_prev_frame_set_file_pos = + trg_data->first_trajectory_frame_set_output_file_pos; + } + else if(trg_data->n_trajectory_frame_sets > trg_data->stride_length + 1) + { + /* FIXME: Currently only working if the previous frame set has its + * long stride pointer already set. This might need some fixing. */ + if(frame_set->long_stride_prev_frame_set_file_pos != -1ULL && + frame_set->long_stride_prev_frame_set_file_pos != 0) + { + trg_init_block(&block); + trg_data->input_file = trg_data->output_file; + + curr_pos = ftell(trg_data->output_file); + fseek(trg_data->output_file, + frame_set->long_stride_prev_frame_set_file_pos, + SEEK_SET); + + if(trg_read_block_header(trg_data, &block) != TRG_SUCCESS) + { + printf("Cannot read frame header. %s: %d\n", + __FILE__, __LINE__); + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + /* Read the next frame set from the previous frame set and one + * long stride step back */ + fseek(trg_data->output_file, block.block_contents_size - 4 * + sizeof(uint64_t), SEEK_CUR); + if(fread(&frame_set->long_stride_prev_frame_set_file_pos, + sizeof(frame_set->long_stride_prev_frame_set_file_pos), + 1, trg_data->output_file) == 0) + { + printf("Cannot read block. %s: %d\n", __FILE__, __LINE__); + trg_destroy_block(&block); + trg_data->input_file = temp; + return(TRG_CRITICAL); + } + + if(trg_data->endianness_64 != TRG_BIG_ENDIAN_64) + { + if(trg_swap_byte_order_64(trg_data, + &frame_set-> + long_stride_prev_frame_set_file_pos) + != TRG_SUCCESS) + { + printf("Cannot swap byte order to get big endian. %s: %d\n", + __FILE__, __LINE__); + } + } + + trg_destroy_block(&block); + trg_data->input_file = temp; + fseek(trg_data->output_file, curr_pos, SEEK_SET); + } + } + + frame_set->first_frame = first_frame; + frame_set->n_frames = n_frames; +// frame_set->n_particle_data_blocks = 0; +// frame_set->n_data_blocks = 0; + + if(trg_data->first_trajectory_frame_set_output_file_pos == -1ULL || + trg_data->first_trajectory_frame_set_output_file_pos == 0) + { + trg_data->first_trajectory_frame_set_output_file_pos = + trg_data->current_trajectory_frame_set_output_file_pos; + } + /* FIXME: Should check the frame number instead of the file_pos, + * in case frame sets are not in order */ + if(trg_data->last_trajectory_frame_set_output_file_pos == -1ULL || + trg_data->last_trajectory_frame_set_output_file_pos == 0 || + trg_data->last_trajectory_frame_set_output_file_pos < + trg_data->current_trajectory_frame_set_output_file_pos) + { + trg_data->last_trajectory_frame_set_output_file_pos = + trg_data->current_trajectory_frame_set_output_file_pos; + } + + return(TRG_SUCCESS); +} + +/* UNTESTED */ +trg_function_status trg_add_data_block(struct trg_trajectory *trg_data, + const uint64_t id, + const char *block_name, + const char datatype, + const uint64_t n_frames, + const uint64_t n_values_per_frame, + const uint64_t stride_length, + uint64_t codec_id, + void *new_data) +{ + int i, j, block_index, size, len; + struct trg_trajectory_frame_set *frame_set; + struct trg_data *data; + void *orig; + trg_block_type block_type_flag; + + frame_set = &trg_data->current_trajectory_frame_set; + + if(trg_data->current_trajectory_frame_set_output_file_pos > 0) + { + block_type_flag = TRG_TRAJECTORY_BLOCK; + } + else + { + block_type_flag = TRG_NON_TRAJECTORY_BLOCK; + } + + block_index = -1; + /* See if there is already a data block of this ID */ + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + for(i = frame_set->n_data_blocks; i-- ;) + { + data = &frame_set->tr_data[i]; + if(data->block_id == id) + { + block_index = i; + break; + } + } + } + else + { + for(i = trg_data->n_data_blocks; i-- ;) + { + data = &trg_data->non_tr_data[i]; + if(data->block_id == id) + { + block_index = i; + break; + } + } + } + + /* Otherwise create a data block */ + if(block_index == -1) + { + if(trg_create_data_block(trg_data, block_type_flag) != + TRG_SUCCESS) + { + printf("Cannot create data block. %s: %d\n", __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + data = &frame_set->tr_data[frame_set->n_data_blocks - 1]; + } + else + { + data = &trg_data->non_tr_data[trg_data->n_data_blocks - 1]; + } + data->block_id = id; + + data->block_name = malloc(strlen(block_name) + 1); + if(!data->block_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + (int)strlen(block_name)+1, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->block_name, block_name, strlen(block_name) + 1); + + data->datatype = datatype; + + data->stride_length = 1; + data->values = 0; + data->n_values_per_frame = 0; + data->n_frames = 0; + data->codec_id = codec_id; + } + + /* Allocate memory */ + if(!data->values || data->n_frames != n_frames || + data->n_values_per_frame != n_values_per_frame) + { + if(trg_allocate_data_mem(trg_data, data, n_frames, + n_values_per_frame) != + TRG_SUCCESS) + { + printf("Cannot allocate data memory. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + switch(datatype) + { + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_CHAR_DATA: + size = sizeof(char); + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + break; + } + + orig = new_data; + + for(i = 0; i < n_frames; i++) + { + for(j = 0; j < n_values_per_frame; j++) + { + /* FIXME: Not optimal to have this switch statement inside + * the for loops. Check if it should be moved outside. */ + switch(datatype) + { + case TRG_CHAR_DATA: + len = min(strlen(new_data) + 1, + TRG_MAX_STR_LEN); + if(data->values[i][j].c) + { + free(data->values[i][j].c); + } + data->values[i][j].c = malloc(len); + if(!data->values[i][j].c) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + len, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->values[i][j].c, + new_data, len); + new_data += len; + break; + case TRG_INT_DATA: + memcpy(&data->values[i][j].i, + new_data, size); + new_data += size; + break; + case TRG_FLOAT_DATA: + memcpy(&data->values[i][j].f, + new_data, size); + new_data += size; + break; + case TRG_DOUBLE_DATA: + default: + memcpy(&data->values[i][j].d, + new_data, size); + new_data += size; + } + } + } + + new_data = orig; + + return(TRG_SUCCESS); +} + + +trg_function_status trg_add_particle_data_block(struct trg_trajectory *trg_data, + const uint64_t id, + const char *block_name, + const char datatype, + const uint64_t n_frames, + const uint64_t n_values_per_frame, + const uint64_t stride_length, + const uint64_t first_particle_number, + const uint64_t n_particles, + uint64_t codec_id, + void *new_data) +{ + int i, j, k, block_index, size, len; + uint64_t tot_n_particles; + struct trg_trajectory_frame_set *frame_set; + struct trg_particle_data *data; + void *orig; + trg_block_type block_type_flag; + + frame_set = &trg_data->current_trajectory_frame_set; + + if(trg_data->current_trajectory_frame_set_output_file_pos > 0) + { + block_type_flag = TRG_TRAJECTORY_BLOCK; + } + else + { + block_type_flag = TRG_NON_TRAJECTORY_BLOCK; + } + + block_index = -1; + /* See if there is already a data block of this ID */ + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + for(i = frame_set->n_particle_data_blocks; i-- ;) + { + data = &frame_set->tr_particle_data[i]; + if(data->block_id == id) + { + block_index = i; + break; + } + } + } + else + { + for(i = trg_data->n_particle_data_blocks; i-- ;) + { + data = &trg_data->non_tr_particle_data[i]; + if(data->block_id == id) + { + block_index = i; + break; + } + } + } + + /* Otherwise create a data block */ + if(block_index == -1) + { + if(trg_create_particle_data_block(trg_data, block_type_flag) != + TRG_SUCCESS) + { + printf("Cannot create particle data block. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + if(block_type_flag == TRG_TRAJECTORY_BLOCK) + { + data = &frame_set->tr_particle_data[frame_set-> + n_particle_data_blocks - 1]; + } + else + { + data = &trg_data->non_tr_particle_data[trg_data-> + n_particle_data_blocks - 1]; + } + data->block_id = id; + + data->block_name = malloc(strlen(block_name) + 1); + if(!data->block_name) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + (int)strlen(block_name)+1, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->block_name, block_name, strlen(block_name) + 1); + + data->datatype = datatype; + + data->stride_length = 1; + data->values = 0; + data->n_values_per_frame = 0; + data->n_frames = 0; + data->codec_id = codec_id; + } + + if(block_type_flag == TRG_TRAJECTORY_BLOCK && trg_data->var_num_atoms_flag) + { + tot_n_particles = frame_set->n_particles; + } + else + { + tot_n_particles = trg_data->n_particles; + } + + /* Allocate memory */ + if(!data->values || data->n_frames != n_frames || + data->n_values_per_frame != n_values_per_frame) + { + if(trg_allocate_particle_data_mem(trg_data, data, n_frames, + tot_n_particles, + n_values_per_frame) != + TRG_SUCCESS) + { + printf("Cannot allocate particle data memory. %s: %d\n", + __FILE__, __LINE__); + return(TRG_CRITICAL); + } + } + + switch(datatype) + { + case TRG_FLOAT_DATA: + size = sizeof(float); + break; + case TRG_INT_DATA: + size = sizeof(uint64_t); + break; + case TRG_CHAR_DATA: + size = sizeof(char); + case TRG_DOUBLE_DATA: + default: + size = sizeof(double); + break; + } + + orig = new_data; + + for(i = 0; i < n_frames; i++) + { + for(j = first_particle_number; j < n_particles; j++) + { + for(k = 0; k < n_values_per_frame; k++) + { + /* FIXME: Not optimal to have this switch statement inside + * the for loops. Check if it should be moved outside. */ + switch(datatype) + { + case TRG_CHAR_DATA: + len = min(strlen(new_data) + 1, + TRG_MAX_STR_LEN); + if(data->values[i][j][k].c) + { + free(data->values[i][j][k].c); + } + data->values[i][j][k].c = malloc(len); + if(!data->values[i][j][k].c) + { + printf("Cannot allocate memory (%d bytes). %s: %d\n", + len, __FILE__, __LINE__); + return(TRG_CRITICAL); + } + strncpy(data->values[i][j][k].c, + new_data, len); + new_data += len; + break; + case TRG_INT_DATA: + memcpy(&data->values[i][j][k].i, + new_data, size); + new_data += size; + break; + case TRG_FLOAT_DATA: + memcpy(&data->values[i][j][k].f, + new_data, size); + new_data += size; + break; + case TRG_DOUBLE_DATA: + default: + memcpy(&data->values[i][j][k].d, + new_data, size); + new_data += size; + } + } + } + } + + new_data = orig; + + return(TRG_SUCCESS); +} + + +trg_function_status trg_read_next_traj_block(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_write_next_traj_block(struct trg_trajectory *trg_data, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_read_traj_block(struct trg_trajectory *trg_data, + uint64_t block_id, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_write_traj_block(struct trg_trajectory *trg_data, + uint64_t block_id, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_read_frame_nr(struct trg_trajectory *trg_data, + uint64_t frame_nr, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_write_frame_nr(struct trg_trajectory *trg_data, + uint64_t frame_nr, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_read_frame_nrs(struct trg_trajectory *trg_data, + uint64_t start_frame_nr, + uint64_t end_frame_nr, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_write_frame_nrs(struct trg_trajectory *trg_data, + uint64_t start_frame_nr, + uint64_t end_frame_nr, + trg_close_file_flag close_file) +{ + /* STUB */ + return(TRG_SUCCESS); +} + +trg_function_status trg_read_frame_set_nr(struct trg_trajectory *trg_data, + uint64_t frame_set_nr) +{ + return(TRG_SUCCESS); +} + + +trg_function_status trg_get_time_str(struct trg_trajectory *trg_data, + char *time) +{ + struct tm *time_data; + time_t secs; + + secs = trg_data->time; + + time_data = localtime(&secs); // Returns a statically allocated variable. + snprintf(time, TRG_MAX_DATE_STR_LEN, + "%4d-%02d-%02d %02d:%02d:%02d", + time_data->tm_year+1900, time_data->tm_mon+1, time_data->tm_mday, + time_data->tm_hour, time_data->tm_min, time_data->tm_sec); + + return(TRG_SUCCESS); +} + diff --git a/src/lib/trg_io.h b/src/lib/trg_io.h new file mode 100644 index 0000000..7bdca5e --- /dev/null +++ b/src/lib/trg_io.h @@ -0,0 +1,760 @@ +#ifndef _TRGIO_H +#define _TRGIO_H 1 + +#include <stdio.h> +#include <inttypes.h> + +#define TRG_VERSION 1 + +#define TRG_UNCOMPRESSED 0ULL +#define TRG_XTC 1ULL +#define TRG_TNG_POSITIONS 2ULL +#define TRG_TNG_VELOCITIES 3ULL +#define TRG_TNG_FORCES 4ULL + +#define TRG_PARTICLE_DEPENDENT 1 +#define TRG_FRAME_DEPENDENT 2 + +// #define TRG_MAX_BLOCK_PARTICLES 1000 +#define TRG_MAX_DATE_STR_LEN 24 +#define TRG_HASH_LEN 16 +#define TRG_MAX_STR_LEN 1024 + +typedef enum {TRG_BIG_ENDIAN_32, + TRG_LITTLE_ENDIAN_32, + TRG_BYTE_PAIR_SWAP_32} trg_endianness_32; + +typedef enum {TRG_BIG_ENDIAN_64, + TRG_LITTLE_ENDIAN_64, + TRG_QUAD_SWAP_64, + TRG_BYTE_PAIR_SWAP_64, + TRG_BYTE_SWAP_64} trg_endianness_64; + + +#define min(a,b) \ + ({ __typeof__ (a) _a = (a); \ + __typeof__ (b) _b = (b); \ + _a < _b ? _a : _b; }) + +#define max(a,b) \ + ({ __typeof__ (a) _a = (a); \ + __typeof__ (b) _b = (b); \ + _a > _b ? _a : _b; }) + +typedef enum {TRG_NON_TRAJECTORY_BLOCK, TRG_TRAJECTORY_BLOCK} trg_block_type; + +typedef enum {TRG_ENDIANNESS_AND_STRING_LENGTH, + TRG_GENERAL_INFO, + TRG_MOLECULES, + TRG_TRAJECTORY_IDS_AND_NAMES, + TRG_TRAJECTORY_FRAME_SET, + TRG_BLOCK_TABLE_OF_CONTENTS, + TRG_PARTICLE_MAPPING} trg_non_trajectory_block_ids; + +typedef enum {TRG_TRAJ_BOX_SHAPE = 10000, + TRG_TRAJ_POSITIONS, + TRG_TRAJ_VELOCITIES, + TRG_TRAJ_FORCES} trg_trajectory_block_ids; + + +typedef enum {TRG_NON_PARTICLE_BLOCK_DATA, + TRG_PARTICLE_BLOCK_DATA} trg_particle_block_data; + +/*typedef enum {TRG_NO_HASH, TRG_OTHER_HASH, + TRG_MD5_HASH, TRG_MD6_HASH, + TRG_SHA0_HASH, TRG_SHA1_HASH, + TRG_SHA256_HASH, TRG_SHA512_HASH} trg_hash_type;*/ + +typedef enum {FALSE, TRUE} trg_bool; + +typedef enum {TRG_KEEP_FILE_OPEN, TRG_CLOSE_FILE} trg_close_file_flag; + +typedef enum {TRG_CONSTANT_N_ATOMS, TRG_VARIABLE_N_ATOMS} trg_variable_n_atoms_flag; + +typedef enum {TRG_SUCCESS, TRG_FAILURE, TRG_CRITICAL} trg_function_status; + +typedef enum {TRG_NORMAL_WRITE, TRG_COPY_EXISTING} write_mode; + +typedef enum {TRG_CHAR_DATA, + TRG_INT_DATA, + TRG_FLOAT_DATA, + TRG_DOUBLE_DATA} trg_data_type; + + + +#ifdef __cplusplus +extern "C" +{ +#endif + +struct trg_bond { + uint64_t from_atom_id; /* One of the atoms of the bond */ + uint64_t to_atom_id; /* The other atom of the bond */ +}; + +struct trg_atom { + struct trg_residue *residue; /* The molecule containing this atom */ + uint64_t id; /* A unique (per molecule) ID number of the atom */ + char *atom_type; /* The atom_type (depending on the forcefield) */ + char *name; /* The name of the atom */ +}; + +struct trg_residue { + struct trg_chain *chain; /* The chain containing this residue */ + uint64_t id; /* A unique (per chain) ID number of the residue */ + char *name; /* The name of the residue */ + uint64_t n_atoms; /* The number of atoms in the residue */ + struct trg_atom *atoms; /* A list of atoms in the residue */ +}; + +struct trg_chain { + struct trg_molecule *molecule; /* The molecule containing this chain */ + uint64_t id; /* A unique (per molecule) ID number of the chain */ + char *name; /* The name of the chain */ + uint64_t n_residues; /* The number of residues in the chain */ + struct trg_residue *residues; /* A list of residues in the chain */ +}; + +struct trg_molecule { + uint64_t id; /* A unique ID number of the molecule */ + uint64_t quaternary_str; /* Quaternary structure of the molecule. + 1 => monomeric + 2 => dimeric + 3 => trimeric + etc */ + uint64_t n_chains; /* The number of chains in the molecule */ + uint64_t n_residues; /* The number of residues in the molecule */ + uint64_t n_atoms; /* The number of atoms in the molecule */ + uint64_t n_bonds; /* The number of bonds in the molecule. + If the bonds are not specified this + value can be 0. */ + char *name; /* The name of the molecule */ + struct trg_chain *chains; /* A list of chains in the molecule */ + struct trg_residue *residues; /* A list of residues in the molecule */ + struct trg_atom *atoms; /* A list of the atoms in the molecule */ + struct trg_bond *bonds; /* A list of the bonds in the molecule */ +}; + +struct trg_gen_block { + uint64_t header_contents_size; /* The size of the block header in bytes */ + uint64_t block_contents_size; /* The size of the block contents in bytes */ + uint64_t id; /* The ID of the block to determine its type */ + char hash[TRG_HASH_LEN]; /* The MD5 hash of the block to verify integrity */ + char *name; /* The name of the block */ + uint64_t block_version; /* The library version used to write the block */ + char *header_contents; /* The full block header contents */ + char *block_contents; /* The full block contents */ +}; + +struct trg_frame_set_toc { + uint64_t n_blocks; /* The number of blocks listed in this + table of contents */ + char **block_names; /* A list of block names */ +}; + +struct trg_particle_mapping { + uint64_t num_first_particle; /* The index number of the first particle + in this mapping block */ + uint64_t n_particles; /* The number of particles list in this + mapping block */ + uint64_t *real_particle_numbers; /* the mapping of index numbers to the + real particle numbers in the + trajectory */ +}; + +struct trg_trajectory_frame_set { + struct trg_frame_set_toc contents; /* The table of contents of this frame set */ + uint64_t n_mapping_blocks; /* The number of different particle + mapping blocks present. */ + struct trg_particle_mapping *mappings; /* The atom mappings of this frame set */ + uint64_t first_frame; /* The first frame of this frame set */ + uint64_t n_frames; /* The number of frames in this frame set */ + uint64_t *molecule_cnt_list; /* A list of the number of each molecule + type - only used when using variable + number of atoms */ + uint64_t n_particles; /* The number of particles/atoms - only + used when using variable number of + atoms */ + uint64_t next_frame_set_file_pos; /* The file position of the next frame set */ + uint64_t prev_frame_set_file_pos; /* The file position of the previous + frame set */ + uint64_t long_stride_next_frame_set_file_pos; /* The file position of the frame + set one stride step ahead */ + uint64_t long_stride_prev_frame_set_file_pos; /* The file position of the frame + set one stride step behind */ + + /* The data blocks in a frame set are trajectory data blocks */ + int n_particle_data_blocks; /* The number of data blocks + of particle dependent data */ + struct trg_particle_data *tr_particle_data; /* A list of data blocks + containing particle + dependent data */ + int n_data_blocks; + struct trg_data *tr_data; +}; + +/* Data can be either double, float, int or a string */ +union data_values { + double d; + float f; + int i; + char *c; +}; + +/* FIXME: Should there be a pointer to a trg_gen_block from each data block? */ +struct trg_particle_data { + uint64_t block_id; /* The block ID of the data block + containing this particle data. + This is used to determine the + kind of data that is stored */ + char *block_name; /* The name of the data block. + This is used to determine the + kind of data that is stored */ + trg_data_type datatype; /* The type of data stored. */ + uint64_t first_frame_with_data; /* The first frame number of the + first data point */ + uint64_t n_frames; /* The number of frames in this + frame set */ + uint64_t n_values_per_frame; /* The number of values stored per + frame */ + uint64_t stride_length; /* The number of frames between + each data point - e.g. when + storing sparse data. */ + uint64_t codec_id; /* ID of the CODEC used for compression + 0 == no compression. */ + double compression_multiplier; /* The multiplier used for getting + integer values for compression */ + union data_values ***values; /* A 3-dimensional array of values, + sized n_frames x n_particles x + n_values_per_frame */ +}; + +struct trg_data { + uint64_t block_id; + char *block_name; /* The name of the data block. + This is used to determine the + kind of data that is stored */ + trg_data_type datatype; /* The type of data stored. */ + uint64_t first_frame_with_data; + uint64_t n_frames; + uint64_t n_values_per_frame; + uint64_t stride_length; + uint64_t codec_id; /* ID of the CODEC used for compression + 0 == no compression. */ + double compression_multiplier; + union data_values **values; /* A 2-dimensional array of values, + sized n_frames x n_values_per_frame */ +}; + + + +struct trg_trajectory { + char *input_file_path; /* The path of the input trajectory file */ + FILE *input_file; /* A handle to the input file */ + long int input_file_pos; /* The reading position of the file */ + long int input_file_len; /* The length of the input file */ + char *output_file_path; /* The path of the output trajectory file */ + FILE *output_file; /* A handle to the output file */ + long int output_file_pos; /* The writing position of the file */ + + trg_endianness_32 endianness_32; /* The endianness of 32 bit values of + the current computer */ + trg_endianness_64 endianness_64; /* The endianness of 64 bit values of + the current computer */ + + char *program_name; /* The name of the program producing + this trajectory */ + char *forcefield_name; /* The forcefield used in the simulations */ + char *user_name; /* The name of the user running the + simulations */ + uint64_t time; /* The time (n seconds since 1970) when + the file was created */ + char *computer_name; /* The name of the computer on which the + simulations were performed */ + char *pgp_signature; /* The PGP signature of the user creating + the file. */ + char var_num_atoms_flag; /* A flag indicating if the number of atoms + can vary throughout the simulation, e.g. + using a grand canonical ensemble */ + uint64_t frame_set_n_frames; /* The number of frames in a frame set. + It is allowed to have frame sets with + fewer frames, but this will help searching + for specific frames */ + uint64_t stride_length; /* The number of frame sets in a long stride + step */ + + uint64_t n_molecules; /* The number of different kinds of + molecules in the trajectory */ + struct trg_molecule *molecules; /* A list of molecules in the trajectory */ + uint64_t *molecule_cnt_list; /* A list of the count of each molecule - + if using variable number of particles + this will be specified in each frame set */ + uint64_t n_particles; /* The total number of particles/atoms. If + using variable number of particles this + will be specified in each frame set */ + + uint64_t n_id_name_pairs; /* The number of ID-name + pairs */ + struct trg_block_id_name_pair *id_name_pairs; /* A list of ID-name + pairs */ + + uint64_t first_trajectory_frame_set_input_file_pos; /* The pos in the src + file of the first + frame set */ + uint64_t first_trajectory_frame_set_output_file_pos; /* The pos in the dest + file of the first + frame set */ + uint64_t last_trajectory_frame_set_input_file_pos; /* The pos in the src + file of the last + frame set */ + uint64_t last_trajectory_frame_set_output_file_pos; /* The pos in the dest + file of the last + frame set */ + struct trg_trajectory_frame_set current_trajectory_frame_set; /* The currently + active frame + set */ + long int current_trajectory_frame_set_input_file_pos; /* The pos in the src + file of the current + frame set */ + long int current_trajectory_frame_set_output_file_pos;/* The pos in the dest + file of the current + frame set */ + uint64_t n_trajectory_frame_sets; /* The number of frame sets + in the trajectory */ + + uint64_t n_trajectory_blocks; /* The number of trajectory blocks + in the file */ + uint64_t n_non_trajectory_blocks; /* The number of non-trajectory blocks + in the file */ + struct trg_gen_block non_trajectory_blocks[32]; /* A list of non-trajectory + blocks */ + + /* These data blocks are non-trajectory data blocks */ + int n_particle_data_blocks; /* The number of non-frame dependent + particle dependent data blocks */ + struct trg_particle_data *non_tr_particle_data; /* A list of data blocks containing particle + dependent data */ + + int n_data_blocks; /* The number of frame and particle independent + data blocks */ + struct trg_data *non_tr_data; /* A list of frame and particle indepdendent data blocks */ +}; + + + +/* Setup a trajectory data container. + *trg_data is a pointer to pre-allocated memory. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_init_trajectory(struct trg_trajectory *trg_data); + +/* Clean up a trajectory data container. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + All allocated memory in the data structure is freed. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_destroy_trajectory(struct trg_trajectory *trg_data); + +/* Set the name of the input file. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_input_file(struct trg_trajectory *trg_data, + const char *file_name); + +/* Set the name of the output file. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_output_file(struct trg_trajectory *trg_data, + const char *file_name); + +/* Set the program name used when creating the trajectory. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + *new_name is a pointer to the string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_program_name(struct trg_trajectory *trg_data, + const char *new_name); + +/* Set the name of the forcefield used in the trajectory. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + *new_name is a pointer to the string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_forcefield_name(struct trg_trajectory *trg_data, + const char *new_name); + +/* Set the name of the user creating the trajectory. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + *new_name is a pointer to the string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_user_name(struct trg_trajectory *trg_data, + const char *new_name); + +/* Set the name of the computer used when creating the trajectory. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + *new_name is a pointer to the string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_computer_name(struct trg_trajectory *trg_data, + const char *new_name); + +/* Set the pgp_signature of the trajectory. + *trg_data is a pointer to pre-allocated memory containing trajectory data. + *signature is a pointer to the string containing the pgp_signature. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_signature(struct trg_trajectory *trg_data, + const char *signature); + +/* Setup a molecule container. + *molecule is a pointer to pre-allocated memory. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_init_molecule(struct trg_molecule *molecule); + +/* Clean up a molecule container. + *molecule is a pointer to pre-allocated memory containing a molecule. + All allocated memory in the data structure is freed. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_destroy_molecule(struct trg_molecule *molecule); + +/* Setup a data block. + *block is a pointer to pre-allocated memory. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_init_block(struct trg_gen_block *block); + +/* Clean up a data block. + *block is a pointer to pre-allocated memory. + All allocated memory in the data structure is freed. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_destroy_block(struct trg_gen_block *block); + +/* Set the name of a data block. + *trg_data is a pointer to the trajectory containing the block. + *block is a pointer to the block to rename. + *new_name is a string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_block_name(struct trg_trajectory *trg_data, + struct trg_gen_block *block, + const char *new_name); + +/* Add a molecule to the trajectory. + *trg_data is a pointer to the trajectory containing the block. + *name is a pointer to the string containing the name of the new molecule. + **molecule is a pointer to a pointer to the newly created molecule. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_add_molecule(struct trg_trajectory *trg_data, + const char *name, + struct trg_molecule **molecule); + +/* Set the name of a molecule. + *trg_data is a pointer to the trajectory containing the molecule. + *molecule is a pointer to the molecule to rename. + *new_name is a string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_molecule_name(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + const char *new_name); + +/* Set the count of a molecule. + *trg_data is a pointer to the trajectory containing the molecule. + *molecule is a pointer to the molecule to rename. + *cnt is a pointer to the variable to be populated with the count. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_get_molecule_cnt(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + uint64_t *cnt); + +/* Set the count of a molecule. + *trg_data is a pointer to the trajectory containing the molecule. + *molecule is a pointer to the molecule to rename. + cnt is the number of instances of this molecule. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_molecule_cnt(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + uint64_t cnt); + +/* Add a chain to a molecule. + *trg_data is a pointer to the trajectory containing the molecule. + *molecule is a pointer to the molecule to add a chain to. + *name is a pointer to a string containing the name of the chain. + **chain is a pointer to a pointer to the newly created chain. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_add_chain_to_molecule(struct trg_trajectory *trg_data, + struct trg_molecule *molecule, + const char *name, + struct trg_chain **chain); + +/* Set the name of a chain. + *trg_data is a pointer to the trajectory containing the atom. + *chain is a pointer to the chain to rename. + *new_name is a string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_chain_name(struct trg_trajectory *trg_data, + struct trg_chain *chain, + const char *new_name); + +/* Add a residue to a chain. + *trg_data is a pointer to the trajectory containing the chain. + *chain is a pointer to the chain to add a residue to. + *name is a pointer to a string containing the name of the residue. + **residue is a pointer to a pointer to the newly created residue. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_add_residue_to_chain(struct trg_trajectory *trg_data, + struct trg_chain *chain, + const char *name, + struct trg_residue **residue); + +/* Set the name of a residue. + *trg_data is a pointer to the trajectory containing the atom. + *residue is a pointer to the residue to rename. + *new_name is a string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_residue_name(struct trg_trajectory *trg_data, + struct trg_residue *residue, + const char *new_name); + +/* Add an atom to a residue. + *trg_data is a pointer to the trajectory containing the residue. + *residue is a pointer to the residue to add an atom to. + *atom_name is a pointer to a string containing the name of the atom. + *atom_type is a pointer to a string containing the atom type of the atom. + **atom is a pointer to a pointer to the newly created atom. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_add_atom_to_residue(struct trg_trajectory *trg_data, + struct trg_residue *residue, + const char *atom_name, + const char *atom_type, + struct trg_atom **atom); + +/* Set the name of an atom. + *trg_data is a pointer to the trajectory containing the atom. + *atom is a pointer to the atom to rename. + *new_name is a string containing the wanted name. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_atom_name(struct trg_trajectory *trg_data, + struct trg_atom *atom, + const char *new_name); + +/* Set the atom type of an atom. + *trg_data is a pointer to the trajectory containing the atom. + *atom is a pointer to the atom to change. + *new_name is a string containing the atom type. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_set_atom_type(struct trg_trajectory *trg_data, + struct trg_atom *atom, + const char *new_type); + +/* Read the header blocks from the input_file of trg_data. + The trajectory blocks must be read separately and iteratively in chunks + to fit in memory. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_file_headers(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + +/* Write the header blocks to the output_file of trg_data. + The trajectory blocks must be written separately and iteratively in chunks + to fit in memory. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (output_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_file_headers(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + + +/* Read one (the next) block (of any kind) from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + *block_data is a pointer to the struct which will be populated with the + data. If block_data->input_file_pos > 0 it is the position from where the reading + starts otherwise it starts from the current position. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor + error has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_next_block(struct trg_trajectory *trg_data, + struct trg_gen_block *block_data, + trg_close_file_flag close_file); + +/* Write one block (of any kind) to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to read from. If the file (output_file) is not open it will be opened. + *block_data is a pointer to the struct containing the + data. If block_data->output_file_pos > 0 it is the position from where the writing + starts otherwise it starts from the current position. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +// trg_function_status trg_write_block(struct trg_trajectory *trg_data, +// struct trg_gen_block *block_data, +// trg_close_file_flag close_file); + + +/* Read one (the next) frame set, including toc, mapping and related data blocks + from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_next_frame_set(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + +/* Write one (the next) frame set, including toc, mapping and related data blocks + to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (output_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_frame_set(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + +trg_function_status trg_new_frame_set(struct trg_trajectory *trg_data, + const uint64_t first_frame, + const uint64_t n_frames); + +trg_function_status trg_add_particle_data_block(struct trg_trajectory *trg_data, + const uint64_t id, + const char *block_name, + const char datatype, + const uint64_t n_frames, + const uint64_t n_values_per_frame, + const uint64_t stride_length, + const uint64_t first_particle_number, + const uint64_t n_particles, + uint64_t codec_id, + void *new_data); + + +/* Read one (the next) trajectory block from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_next_traj_block(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + +/* Write one (the next) trajectory block to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (output_file) is not open it will be opened. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_next_traj_block(struct trg_trajectory *trg_data, + trg_close_file_flag close_file); + +/* Read one trajectory block from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + block_id is the ID of the block to read. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_traj_block(struct trg_trajectory *trg_data, + uint64_t block_id, + trg_close_file_flag close_file); + +/* Write one trajectory block to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (output_file) is not open it will be opened. + block_id is the ID of the block to write. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_traj_block(struct trg_trajectory *trg_data, + uint64_t block_id, + trg_close_file_flag close_file); + +/* Read a requested frame set. + *trg_data is a pointer to trajectory data. trg_data->current_trajectory_frame_set + will be the read frame set. + frame_set_nr is the number of the frame set to return (starting from 0). + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_frame_set_nr(struct trg_trajectory *trg_data, + uint64_t frame_set_nr); + +/* Read one trajectory frame from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + frame_nr is the index number of the frame to read. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_frame_nr(struct trg_trajectory *trg_data, + uint64_t frame_nr, + trg_close_file_flag close_file); + +/* Write one trajectory frame to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (output_file) is not open it will be opened. + frame_nr is the index number of the frame to write. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_frame_nr(struct trg_trajectory *trg_data, + uint64_t frame_nr, + trg_close_file_flag close_file); + +/* Read a number of consecutive trajectory frames from the input_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->input_file_path specifies which + file to read from. If the file (input_file) is not open it will be opened. + start_frame_nr is the index number of the first frame to read. + end_frame_nr is the index number of the last frame to read. + If close_file == TRG_CLOSE_FILE (1) the input_file will be closed after reading the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_read_frame_nrs(struct trg_trajectory *trg_data, + uint64_t start_frame_nr, + uint64_t end_frame_nr, + trg_close_file_flag close_file); + +/* Write a number of consecutive trajectory frames to the output_file of trg_data. + *trg_data is a pointer to trajectory data. trg_data->output_file_path specifies which + file to write to. If the file (input_file) is not open it will be opened. + start_frame_nr is the index number of the first frame to write. + end_frame_nr is the index number of the last frame to write. + If close_file == TRG_CLOSE_FILE (1) the output_file will be closed after writing the data. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_write_frame_nrs(struct trg_trajectory *trg_data, + uint64_t start_frame_nr, + uint64_t end_frame_nr, + trg_close_file_flag close_file); + +/* Get the date and time of initial file creation in ISO format (string). + *trg_data is a pointer to trajectory data. + *time is a pointer to the string in which the date will be stored. Memory + must be reserved beforehand. + Returns TRG_SUCCESS (0) if successful, TRG_FAILURE (1) if a minor error + has occurred or TRG_CRITICAL (2) if a major error has occured. */ +trg_function_status trg_get_time_str(struct trg_trajectory *trg_data, char *time); + +#ifdef __cplusplus +} /* end extern "C" */ +#endif + +#endif /* _TRGIO_H */
\ No newline at end of file diff --git a/src/tests/CMakeLists.txt b/src/tests/CMakeLists.txt new file mode 100644 index 0000000..ea5a644 --- /dev/null +++ b/src/tests/CMakeLists.txt @@ -0,0 +1,5 @@ +include_directories(${TRAJECTORY_SOURCE_DIR}/src/lib) +link_directories(${TRAJECTORY_BINARY_DIR}/src/lib) + +add_executable(trg_testing trg_io_testing.c) +target_link_libraries(trg_testing trg_io)
\ No newline at end of file diff --git a/src/tests/trg_io_testing.c b/src/tests/trg_io_testing.c new file mode 100644 index 0000000..caa74dd --- /dev/null +++ b/src/tests/trg_io_testing.c @@ -0,0 +1,341 @@ +#include <inttypes.h> +#include <stdlib.h> +#include <string.h> +#include "trg_io.h" + + + + +static trg_function_status trg_setup_test_molecules(struct trg_trajectory *traj) +{ + struct trg_molecule *molecule; + struct trg_chain *chain; + struct trg_residue *residue; + struct trg_atom *atom; + uint64_t cnt; +// int i; + + trg_add_molecule(traj, "water", &molecule); + trg_add_chain_to_molecule(traj, molecule, "W", &chain); + trg_add_residue_to_chain(traj, chain, "WAT", &residue); + if(trg_add_atom_to_residue(traj, residue, "O", "O", &atom) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + if(trg_add_atom_to_residue(traj, residue, "HO1", "H", &atom) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + if(trg_add_atom_to_residue(traj, residue, "HO2", "H", &atom) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + trg_set_molecule_cnt(traj, molecule, 200); + trg_get_molecule_cnt(traj, molecule, &cnt); + printf("Created %"PRId64" %s molecules.\n", cnt, molecule->name); +// traj->molecule_cnt_list[traj->n_molecules-1] = 5; +// trg_set_molecule_name(traj, &traj->molecules[1], "ligand"); +// trg_set_molecule_name(traj, &traj->molecules[2], "water"); +// trg_set_molecule_name(traj, &traj->molecules[3], "dummy"); +// traj->molecules[0].id = 0; +// traj->molecules[1].id = 1; +// traj->molecules[2].id = 2; +// traj->molecules[3].id = 3; + +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom1", "type1") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom2", "type1") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom3", "type1") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom4", "type2") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom5", "type2") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom6", "type2") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[0], "atom7", "type3") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[1], "C1", "C") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[1], "O1", "O") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[1], "H11", "H") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[1], "H12", "H") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } +// if(trg_add_atom_to_molecule(traj, &traj->molecules[1], "H13", "H") == TRG_CRITICAL) +// { +// return(TRG_CRITICAL); +// } + + return(TRG_SUCCESS); +} + +static trg_function_status trg_test_read_and_write_file + (struct trg_trajectory *traj) +{ + trg_function_status stat; + + stat = trg_read_file_headers(traj, TRG_KEEP_FILE_OPEN); + if(stat == TRG_CRITICAL) + { + return(stat); + } + stat = trg_write_file_headers(traj, TRG_KEEP_FILE_OPEN); + if(stat == TRG_CRITICAL) + { + return(stat); + } + + while(stat != TRG_CRITICAL && traj->input_file_pos < traj->input_file_len && + traj->current_trajectory_frame_set.next_frame_set_file_pos != -1UL) + { + stat = trg_read_next_frame_set(traj, TRG_KEEP_FILE_OPEN); + if(stat == TRG_CRITICAL) + { + return(stat); + } + stat = trg_write_frame_set(traj, TRG_KEEP_FILE_OPEN); + } + + return(stat); +} + +static trg_function_status trg_test_write_and_read_traj(struct trg_trajectory *traj) +{ + int i, j, k, nr, tot_n_mols, cnt; + float *data, *molpos; + trg_function_status stat; + + /* Create molecules */ + if(trg_setup_test_molecules(traj) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + + if(trg_init_block(&traj->non_trajectory_blocks[traj->n_non_trajectory_blocks]) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + traj->non_trajectory_blocks[traj->n_non_trajectory_blocks].id = TRG_MOLECULES; + if(trg_set_block_name(traj, + &traj->non_trajectory_blocks[traj->n_non_trajectory_blocks++], + "MOLECULES") == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + + if(trg_write_file_headers(traj, TRG_KEEP_FILE_OPEN) == TRG_CRITICAL) + { + return(TRG_CRITICAL); + } + + data = malloc(sizeof(float) * traj->n_particles * + traj->frame_set_n_frames * 3); + if(!data) + { + printf("Cannot allocate memory. %s: %d\n", __FILE__, __LINE__); + return(TRG_CRITICAL); + } + + tot_n_mols = 0; + for(i = 0; i < traj->n_molecules; i++) + { + tot_n_mols += traj->molecule_cnt_list[i]; + } + molpos = malloc(sizeof(float) * tot_n_mols * 3); + + /* Set initial coordinates */ + for(i = 0; i < tot_n_mols; i++) + { + nr = i * 3; + /* Somewhat random coordinates (between 0 and 100), + * but not specifying a random seed */ + molpos[nr] = 100.0 * rand() / (RAND_MAX + 1.0); + molpos[nr+1] = 100.0 * rand() / (RAND_MAX + 1.0); + molpos[nr+2] = 100.0 * rand() / (RAND_MAX + 1.0); + } + + /* Generate 200 frame sets - each with 100 frames (by default) */ + for(i = 0; i < 200; i++) + { + cnt = 0; + for(j = 0; j < traj->frame_set_n_frames; j++) + { + for(k = 0; k < tot_n_mols; k++) + { + nr = k * 3; + /* Move -1 to 1 */ + molpos[nr] += 2 * (rand() / (RAND_MAX + 1.0)) - 1; + molpos[nr+1] += 2 * (rand() / (RAND_MAX + 1.0)) - 1; + molpos[nr+2] += 2 * (rand() / (RAND_MAX + 1.0)) - 1; + + data[cnt++] = molpos[nr]; + data[cnt++] = molpos[nr + 1]; + data[cnt++] = molpos[nr + 2]; + data[cnt++] = molpos[nr] + 1; + data[cnt++] = molpos[nr + 1] + 1; + data[cnt++] = molpos[nr + 2] + 1; + data[cnt++] = molpos[nr] - 1; + data[cnt++] = molpos[nr + 1] - 1; + data[cnt++] = molpos[nr + 2] - 1; + } + } + if(trg_new_frame_set(traj, i * traj->frame_set_n_frames, + traj->frame_set_n_frames) != TRG_SUCCESS) + { + printf("Error creating frame set %d. %s: %d\n", + i, __FILE__, __LINE__); + free(molpos); + free(data); + return(TRG_CRITICAL); + } + + if(trg_add_particle_data_block(traj, TRG_TRAJ_POSITIONS, + "POSITIONS", + TRG_FLOAT_DATA, + traj->frame_set_n_frames, 3, + 1, 0, traj->n_particles, + TRG_UNCOMPRESSED, + data) != TRG_SUCCESS) + { + printf("Error adding data. %s: %d\n", __FILE__, __LINE__); + free(molpos); + free(data); + return(TRG_CRITICAL); + } + if(trg_write_frame_set(traj, TRG_KEEP_FILE_OPEN) != TRG_SUCCESS) + { + printf("Error writing frame set. %s: %d\n", __FILE__, __LINE__); + free(molpos); + free(data); + return(TRG_CRITICAL); + } + } + +// trg_add_ids_names_pair(traj, TRG_TRAJ_VELOCITIES, "BOX SHAPE"); +// trg_add_ids_names_pair(traj, TRG_TRAJ_POSITIONS, "TRAJECTORY POSITIONS"); +// trg_add_ids_names_pair(traj, TRG_TRAJ_VELOCITIES, "TRAJECTORY VELOCITIES"); +// trg_add_ids_names_pair(traj, TRG_TRAJ_FORCES, "TRAJECTORY FORCES"); +// trg_add_ids_names_pair(traj, 11000, "TEST DATA"); + + free(molpos); + free(data); + + trg_destroy_trajectory(traj); + trg_set_input_file(traj, "/tmp/trg_test.trg"); + + stat = trg_read_file_headers(traj, TRG_KEEP_FILE_OPEN); + + while(stat != TRG_CRITICAL && traj->input_file_pos < traj->input_file_len && + traj->current_trajectory_frame_set.next_frame_set_file_pos != -1ULL) + { + stat = trg_read_next_frame_set(traj, TRG_KEEP_FILE_OPEN); + if(stat == TRG_CRITICAL) + { + return(stat); + } + } + + return(stat); +} + +int main() +{ + struct trg_trajectory traj; + char time_str[TRG_MAX_DATE_STR_LEN]; + + if(trg_init_trajectory(&traj) != TRG_SUCCESS) + { + trg_destroy_trajectory(&traj); + printf("Test Init trajectory:\t\t\t\tFailed. %s: %d.\n", + __FILE__, __LINE__); + exit(1); + } + printf("Test Init trajectory:\t\t\t\tSucceeded.\n"); + + trg_get_time_str(&traj, time_str); + + printf("Creation time: %s\n", time_str); + + trg_set_input_file(&traj, "trg_example.trg"); + trg_set_output_file(&traj, "/tmp/trg_example_test.trg"); + +// if(trg_test_endianness(&traj) != TRG_SUCCESS) +// { +// printf("Test failed: Endianness. %s: %d\n", __FILE__, __LINE__); +// } + + if(trg_test_read_and_write_file(&traj) == TRG_CRITICAL) + { + printf("Test Read and write file:\t\t\tFailed. %s: %d\n", + __FILE__, __LINE__); + } + else + { + printf("Test Read and write file:\t\t\tSucceeded.\n"); + } + + if(trg_destroy_trajectory(&traj) == TRG_CRITICAL || + trg_init_trajectory(&traj) == TRG_CRITICAL) + { + printf("Test Destroy and init trajectory:\t\tFailed. %s: %d\n", + __FILE__, __LINE__); + } + else + { + printf("Test Destroy and init trajectory:\t\tSucceeded.\n"); + } + + + trg_set_output_file(&traj, "/tmp/trg_test.trg"); + + if(trg_test_write_and_read_traj(&traj) == TRG_CRITICAL) + { + printf("Test Write and read file:\t\t\tFailed. %s: %d\n", + __FILE__, __LINE__); + } + else + { + printf("Test Write and read file:\t\t\tSucceeded.\n"); + } + + if(trg_destroy_trajectory(&traj) == TRG_CRITICAL) + { + printf("Test Destroy trajectory:\t\t\tFailed. %s: %d.\n", + __FILE__, __LINE__); + exit(1); + } + else + { + printf("Test Destroy trajectory:\t\t\tSucceeded.\n"); + } + + + printf("Tests finished\n"); + + exit(0); +}
\ No newline at end of file |