Bug Summary

File:external/tng_io/src/compression/tng_compress.c
Location:line 383, column 5
Description:Null pointer passed as an argument to a 'nonnull' parameter

Annotated Source Code

1/* This code is part of the tng compression routines.
2 *
3 * Written by Daniel Spangberg
4 * Copyright (c) 2010, 2013, The GROMACS development team.
5 *
6 *
7 * This program is free software; you can redistribute it and/or
8 * modify it under the terms of the Revised BSD License.
9 */
10
11#include <stdlib.h>
12#include <string.h>
13#include <math.h>
14#include "../../include/compression/tng_compress.h"
15#include "../../include/compression/coder.h"
16#include "../../include/compression/fixpoint.h"
17
18/* Please see tng_compress.h for info on how to call these routines. */
19
20/* This becomes TNGP for positions (little endian) and TNGV for velocities. In ASCII. */
21#define MAGIC_INT_POS0x50474E54 0x50474E54
22#define MAGIC_INT_VEL0x56474E54 0x56474E54
23
24#define SPEED_DEFAULT2 2 /* Default to relatively fast compression. For very good compression it makes sense to
25 choose speed=4 or speed=5 */
26
27#define PRECISION(hi,lo)(Ptngc_i32x2_to_d(hi,lo)) (Ptngc_i32x2_to_d(hi,lo))
28
29#define MAX_FVAL2147483647. 2147483647.
30
31static int verify_input_data(double *x, int natoms, int nframes, double precision)
32{
33 int iframe, i, j;
34 for (iframe=0; iframe<nframes; iframe++)
35 for (i=0; i<natoms; i++)
36 for (j=0; j<3; j++)
37 if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL2147483647.)
38 goto error;
39 return 0;
40 error:
41#if 0
42 for (iframe=0; iframe<nframes; iframe++)
43 for (i=0; i<natoms; i++)
44 for (j=0; j<3; j++)
45 if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL2147483647.)
46 printf("ERROR. Too large value: %d %d %d: %g %g %g\n",iframe,i,j,x[iframe*natoms*3+i*3+j],precision,x[iframe*natoms*3+i*3+j]/precision/MAX_FVAL2147483647.);
47#endif
48 return 1;
49}
50
51static int verify_input_data_float(float *x, int natoms, int nframes, float precision)
52{
53 int iframe, i, j;
54 for (iframe=0; iframe<nframes; iframe++)
55 for (i=0; i<natoms; i++)
56 for (j=0; j<3; j++)
57 if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL2147483647.)
58 goto error;
59 return 0;
60 error:
61#if 0
62 for (iframe=0; iframe<nframes; iframe++)
63 for (i=0; i<natoms; i++)
64 for (j=0; j<3; j++)
65 if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL2147483647.)
66 printf("ERROR. Too large value: %d %d %d: %g %g %g\n",iframe,i,j,x[iframe*natoms*3+i*3+j],precision,x[iframe*natoms*3+i*3+j]/precision/MAX_FVAL2147483647.);
67#endif
68 return 1;
69}
70
71static int quantize(double *x, int natoms, int nframes,
72 double precision,
73 int *quant)
74{
75 int iframe, i, j;
76 for (iframe=0; iframe<nframes; iframe++)
77 for (i=0; i<natoms; i++)
78 for (j=0; j<3; j++)
79 quant[iframe*natoms*3+i*3+j]=(int)floor((x[iframe*natoms*3+i*3+j]/precision)+0.5);
80 return verify_input_data(x,natoms,nframes,precision);
81}
82
83static int quantize_float(float *x, int natoms, int nframes,
84 float precision,
85 int *quant)
86{
87 int iframe, i, j;
88 for (iframe=0; iframe<nframes; iframe++)
89 for (i=0; i<natoms; i++)
90 for (j=0; j<3; j++)
91 quant[iframe*natoms*3+i*3+j]=(int)floor((x[iframe*natoms*3+i*3+j]/precision)+0.5);
92 return verify_input_data_float(x,natoms,nframes,precision);
93}
94
95static void quant_inter_differences(int *quant, int natoms, int nframes,
96 int *quant_inter)
97{
98 int iframe, i, j;
99 /* The first frame is used for absolute positions. */
100 for (i=0; i<natoms; i++)
101 for (j=0; j<3; j++)
102 quant_inter[i*3+j]=quant[i*3+j];
103 /* For all other frames, the difference to the previous frame is used. */
104 for (iframe=1; iframe<nframes; iframe++)
105 for (i=0; i<natoms; i++)
106 for (j=0; j<3; j++)
107 quant_inter[iframe*natoms*3+i*3+j]=quant[iframe*natoms*3+i*3+j]-quant[(iframe-1)*natoms*3+i*3+j];
108}
109
110static void quant_intra_differences(int *quant, int natoms, int nframes,
111 int *quant_intra)
112{
113 int iframe, i, j;
114 for (iframe=0; iframe<nframes; iframe++)
115 {
116 /* The first atom is used with its absolute position. */
117 for (j=0; j<3; j++)
118 quant_intra[iframe*natoms*3+j]=quant[iframe*natoms*3+j];
119 /* For all other atoms the intraframe differences are computed. */
120 for (i=1; i<natoms; i++)
121 for (j=0; j<3; j++)
122 quant_intra[iframe*natoms*3+i*3+j]=quant[iframe*natoms*3+i*3+j]-quant[iframe*natoms*3+(i-1)*3+j];
123 }
124}
125
126static void unquantize(double *x, int natoms, int nframes,
127 double precision,
128 int *quant)
129{
130 int iframe, i, j;
131 for (iframe=0; iframe<nframes; iframe++)
132 for (i=0; i<natoms; i++)
133 for (j=0; j<3; j++)
134 x[iframe*natoms*3+i*3+j]=(double)quant[iframe*natoms*3+i*3+j]*precision;
135}
136
137static void unquantize_float(float *x, int natoms, int nframes,
138 float precision,
139 int *quant)
140{
141 int iframe, i, j;
142 for (iframe=0; iframe<nframes; iframe++)
143 for (i=0; i<natoms; i++)
144 for (j=0; j<3; j++)
145 x[iframe*natoms*3+i*3+j]=(float)quant[iframe*natoms*3+i*3+j]*precision;
146}
147
148static void unquantize_inter_differences(double *x, int natoms, int nframes,
149 double precision,
150 int *quant)
151{
152 int iframe, i, j;
153 for (i=0; i<natoms; i++)
154 for (j=0; j<3; j++)
155 {
156 int q=quant[i*3+j]; /* First value. */
157 x[i*3+j]=(double)q*precision;
158 for (iframe=1; iframe<nframes; iframe++)
159 {
160 q+=quant[iframe*natoms*3+i*3+j];
161 x[iframe*natoms*3+i*3+j]=(double)q*precision;
162 }
163 }
164}
165
166static void unquantize_inter_differences_float(float *x, int natoms, int nframes,
167 float precision,
168 int *quant)
169{
170 int iframe, i, j;
171 for (i=0; i<natoms; i++)
172 for (j=0; j<3; j++)
173 {
174 int q=quant[i*3+j]; /* First value. */
175 x[i*3+j]=(float)q*precision;
176 for (iframe=1; iframe<nframes; iframe++)
177 {
178 q+=quant[iframe*natoms*3+i*3+j];
179 x[iframe*natoms*3+i*3+j]=(float)q*precision;
180 }
181 }
182}
183
184static void unquantize_inter_differences_int(int *x, int natoms, int nframes,
185 int *quant)
186{
187 int iframe, i, j;
188 for (i=0; i<natoms; i++)
189 for (j=0; j<3; j++)
190 {
191 int q=quant[i*3+j]; /* First value. */
192 x[i*3+j]=q;
193 for (iframe=1; iframe<nframes; iframe++)
194 {
195 q+=quant[iframe*natoms*3+i*3+j];
196 x[iframe*natoms*3+i*3+j]=q;
197 }
198 }
199}
200
201/* In frame update required for the initial frame if intra-frame
202 compression was used. */
203static void unquant_intra_differences_first_frame(int *quant, int natoms)
204{
205 int i,j;
206 for (j=0; j<3; j++)
207 {
208 int q=quant[j];
209 for (i=1; i<natoms; i++)
210 {
211 q+=quant[i*3+j];
212 quant[i*3+j]=q;
213 }
214 }
215#if 0
216 for (j=0; j<3; j++)
217 for (i=0; i<natoms; i++)
218 {
219 printf("UQ: %d %d %d: %d\n",0,j,i,quant[i*3+j]);
220 }
221#endif
222}
223
224static void unquantize_intra_differences(double *x, int natoms, int nframes,
225 double precision,
226 int *quant)
227{
228 int iframe, i, j;
229#if 0
230 printf("UQ precision=%g\n",precision);
231#endif
232 for (iframe=0; iframe<nframes; iframe++)
233 for (j=0; j<3; j++)
234 {
235 int q=quant[iframe*natoms*3+j];
236 x[iframe*natoms*3+j]=(double)q*precision;
237 for (i=1; i<natoms; i++)
238 {
239 q+=quant[iframe*natoms*3+i*3+j];
240 x[iframe*natoms*3+i*3+j]=(double)q*precision;
241 }
242 }
243}
244
245static void unquantize_intra_differences_float(float *x, int natoms, int nframes,
246 float precision,
247 int *quant)
248{
249 int iframe, i, j;
250 for (iframe=0; iframe<nframes; iframe++)
251 for (j=0; j<3; j++)
252 {
253 int q=quant[iframe*natoms*3+j];
254 x[iframe*natoms*3+j]=(float)q*precision;
255 for (i=1; i<natoms; i++)
256 {
257 q+=quant[iframe*natoms*3+i*3+j];
258 x[iframe*natoms*3+i*3+j]=(float)q*precision;
259 }
260 }
261}
262
263static void unquantize_intra_differences_int(int *x, int natoms, int nframes,
264 int *quant)
265{
266 int iframe, i, j;
267 for (iframe=0; iframe<nframes; iframe++)
268 for (j=0; j<3; j++)
269 {
270 int q=quant[iframe*natoms*3+j];
271 x[iframe*natoms*3+j]=q;
272 for (i=1; i<natoms; i++)
273 {
274 q+=quant[iframe*natoms*3+i*3+j];
275 x[iframe*natoms*3+i*3+j]=q;
276 }
277 }
278}
279
280/* Buffer num 8 bit bytes into buffer location buf */
281static void bufferfix(unsigned char *buf, fix_t v, int num)
282{
283 /* Store in little endian format. */
284 unsigned char c; /* at least 8 bits. */
285 c=(unsigned char)(v & 0xFFU);
286 while (num--)
287 {
288 *buf++=c;
289 v >>= 8;
290 c=(unsigned char)(v & 0xFFU);
291 }
292}
293
294static fix_t readbufferfix(unsigned char *buf, int num)
295{
296 unsigned char b;
297 int shift=0;
298 fix_t f=0UL;
299 int cnt=0;
300 do
301 {
302 b=buf[cnt++];
303 f |= ((fix_t)b & 0xFF)<<shift;
304 shift+=8;
305 } while (--num);
306 return f;
307}
308
309/* Perform position compression from the quantized data. */
310static void compress_quantized_pos(int *quant, int *quant_inter, int *quant_intra,
311 int natoms, int nframes,
312 int speed,
313 int initial_coding, int initial_coding_parameter,
314 int coding, int coding_parameter,
315 fix_t prec_hi, fix_t prec_lo,
316 int *nitems,
317 char *data)
318{
319 int bufloc=0;
320 char *datablock=NULL((void*)0);
1
'datablock' initialized to a null pointer value
→
321 int length=0;
322 /* Information needed for decompression. */
323 if (data)
2
←
Assuming 'data' is non-null
→
3
←
Taking true branch
→
324 bufferfix((unsigned char*)data+bufloc,(fix_t)MAGIC_INT_POS0x50474E54,4);
325 bufloc+=4;
326 /* Number of atoms. */
327 if (data)
4
←
Taking true branch
→
328 bufferfix((unsigned char*)data+bufloc,(fix_t)natoms,4);
329 bufloc+=4;
330 /* Number of frames. */
331 if (data)
5
←
Taking true branch
→
332 bufferfix((unsigned char*)data+bufloc,(fix_t)nframes,4);
333 bufloc+=4;
334 /* Initial coding. */
335 if (data)
6
←
Taking true branch
→
336 bufferfix((unsigned char*)data+bufloc,(fix_t)initial_coding,4);
337 bufloc+=4;
338 /* Initial coding parameter. */
339 if (data)
7
←
Taking true branch
→
340 bufferfix((unsigned char*)data+bufloc,(fix_t)initial_coding_parameter,4);
341 bufloc+=4;
342 /* Coding. */
343 if (data)
8
←
Taking true branch
→
344 bufferfix((unsigned char*)data+bufloc,(fix_t)coding,4);
345 bufloc+=4;
346 /* Coding parameter. */
347 if (data)
9
←
Taking true branch
→
348 bufferfix((unsigned char*)data+bufloc,(fix_t)coding_parameter,4);
349 bufloc+=4;
350 /* Precision. */
351 if (data)
10
←
Taking true branch
→
352 bufferfix((unsigned char*)data+bufloc,prec_lo,4);
353 bufloc+=4;
354 if (data)
11
←
Taking true branch
→
355 bufferfix((unsigned char*)data+bufloc,prec_hi,4);
356 bufloc+=4;
357 /* The initial frame */
358 if ((initial_coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
12
←
Assuming 'initial_coding' is not equal to 5
→
15
←
Taking false branch
→
359 (initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7) ||
13
←
Assuming 'initial_coding' is not equal to 7
→
360 (initial_coding==TNG_COMPRESS_ALGO_POS_XTC310))
14
←
Assuming 'initial_coding' is not equal to 10
→
361 {
362 struct coder *coder=Ptngc_coder_init();
363 length=natoms*3;
364 datablock=(char*)Ptngc_pack_array(coder,quant,&length,
365 initial_coding,initial_coding_parameter,natoms,speed);
366 Ptngc_coder_deinit(coder);
367 }
368 else if ((initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3) ||
16
←
Assuming 'initial_coding' is not equal to 3
→
18
←
Taking false branch
→
369 (initial_coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
17
←
Assuming 'initial_coding' is not equal to 9
→
370 {
371 struct coder *coder=Ptngc_coder_init();
372 length=natoms*3;
373 datablock=(char*)Ptngc_pack_array(coder,quant_intra,&length,
374 initial_coding,initial_coding_parameter,natoms,speed);
375 Ptngc_coder_deinit(coder);
376 }
377 /* Block length. */
378 if (data)
19
←
Taking true branch
→
379 bufferfix((unsigned char*)data+bufloc,(fix_t)length,4);
380 bufloc+=4;
381 /* The actual data block. */
382 if (data)
20
←
Taking true branch
→
383 memcpy(data+bufloc,datablock,length);
21
←
Null pointer passed as an argument to a 'nonnull' parameter
384 free(datablock);
385 bufloc+=length;
386 /* The remaining frames */
387 if (nframes>1)
388 {
389 datablock=NULL((void*)0);
390 /* Inter-frame compression? */
391 if ((coding==TNG_COMPRESS_ALGO_POS_STOPBIT_INTER1) ||
392 (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTER2) ||
393 (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTER8))
394 {
395 struct coder *coder=Ptngc_coder_init();
396 length=natoms*3*(nframes-1);
397 datablock=(char*)Ptngc_pack_array(coder,quant_inter+natoms*3,&length,
398 coding,coding_parameter,natoms,speed);
399 Ptngc_coder_deinit(coder);
400 }
401 /* One-to-one compression? */
402 else if ((coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
403 (coding==TNG_COMPRESS_ALGO_POS_XTC310) ||
404 (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7))
405 {
406 struct coder *coder=Ptngc_coder_init();
407 length=natoms*3*(nframes-1);
408 datablock=(char*)Ptngc_pack_array(coder,quant+natoms*3,&length,
409 coding,coding_parameter,natoms,speed);
410 Ptngc_coder_deinit(coder);
411 }
412 /* Intra-frame compression? */
413 else if ((coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3) ||
414 (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
415 {
416 struct coder *coder=Ptngc_coder_init();
417 length=natoms*3*(nframes-1);
418 datablock=(char*)Ptngc_pack_array(coder,quant_intra+natoms*3,&length,
419 coding,coding_parameter,natoms,speed);
420 Ptngc_coder_deinit(coder);
421 }
422 /* Block length. */
423 if (data)
424 bufferfix((unsigned char*)data+bufloc,(fix_t)length,4);
425 bufloc+=4;
426 if (datablock)
427 {
428 if (data)
429 memcpy(data+bufloc,datablock,length);
430 free(datablock);
431 }
432 bufloc+=length;
433 }
434 *nitems=bufloc;
435}
436
437/* Perform velocity compression from vel into the data block */
438static void compress_quantized_vel(int *quant, int *quant_inter,
439 int natoms, int nframes,
440 int speed,
441 int initial_coding, int initial_coding_parameter,
442 int coding, int coding_parameter,
443 fix_t prec_hi, fix_t prec_lo,
444 int *nitems,
445 char *data)
446{
447 int bufloc=0;
448 char *datablock=NULL((void*)0);
449 int length;
450 /* Information needed for decompression. */
451 if (data)
452 bufferfix((unsigned char*)data+bufloc,(fix_t)MAGIC_INT_VEL0x56474E54,4);
453 bufloc+=4;
454 /* Number of atoms. */
455 if (data)
456 bufferfix((unsigned char*)data+bufloc,(fix_t)natoms,4);
457 bufloc+=4;
458 /* Number of frames. */
459 if (data)
460 bufferfix((unsigned char*)data+bufloc,(fix_t)nframes,4);
461 bufloc+=4;
462 /* Initial coding. */
463 if (data)
464 bufferfix((unsigned char*)data+bufloc,(fix_t)initial_coding,4);
465 bufloc+=4;
466 /* Initial coding parameter. */
467 if (data)
468 bufferfix((unsigned char*)data+bufloc,(fix_t)initial_coding_parameter,4);
469 bufloc+=4;
470 /* Coding. */
471 if (data)
472 bufferfix((unsigned char*)data+bufloc,(fix_t)coding,4);
473 bufloc+=4;
474 /* Coding parameter. */
475 if (data)
476 bufferfix((unsigned char*)data+bufloc,(fix_t)coding_parameter,4);
477 bufloc+=4;
478 /* Precision. */
479 if (data)
480 bufferfix((unsigned char*)data+bufloc,prec_lo,4);
481 bufloc+=4;
482 if (data)
483 bufferfix((unsigned char*)data+bufloc,prec_hi,4);
484 bufloc+=4;
485
486 length=natoms*3;
487 /* The initial frame */
488 if ((initial_coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1) ||
489 (initial_coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3) ||
490 (initial_coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9))
491 {
492 struct coder *coder=Ptngc_coder_init();
493 datablock=(char*)Ptngc_pack_array(coder,quant,&length,
494 initial_coding,initial_coding_parameter,natoms,speed);
495 Ptngc_coder_deinit(coder);
496 }
497 /* Block length. */
498 if (data)
499 bufferfix((unsigned char*)data+bufloc,(fix_t)length,4);
500 bufloc+=4;
501 /* The actual data block. */
502 if (data && datablock)
503 {
504 memcpy(data+bufloc,datablock,length);
505 free(datablock);
506 bufloc+=length;
507 }
508 /* The remaining frames */
509 if (nframes>1)
510 {
511 datablock=NULL((void*)0);
512 /* Inter-frame compression? */
513 if ((coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_INTER2) ||
514 (coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER6) ||
515 (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_INTER8))
516 {
517 struct coder *coder=Ptngc_coder_init();
518 length=natoms*3*(nframes-1);
519 datablock=(char*)Ptngc_pack_array(coder,quant_inter+natoms*3,&length,
520 coding,coding_parameter,natoms,speed);
521 Ptngc_coder_deinit(coder);
522 }
523 /* One-to-one compression? */
524 else if ((coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1) ||
525 (coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3) ||
526 (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9))
527 {
528 struct coder *coder=Ptngc_coder_init();
529 length=natoms*3*(nframes-1);
530 datablock=(char*)Ptngc_pack_array(coder,quant+natoms*3,&length,
531 coding,coding_parameter,natoms,speed);
532 Ptngc_coder_deinit(coder);
533 }
534 /* Block length. */
535 if (data)
536 bufferfix((unsigned char*)data+bufloc,(fix_t)length,4);
537 bufloc+=4;
538 if (data)
539 memcpy(data+bufloc,datablock,length);
540 free(datablock);
541 bufloc+=length;
542 }
543 *nitems=bufloc;
544}
545
546static int determine_best_coding_stop_bits(struct coder *coder,int *input, int *length,
547 int *coding_parameter, int natoms)
548{
549 int bits;
550 unsigned char *packed;
551 int best_length=0;
552 int new_parameter=-1;
553 int io_length;
554 for (bits=1; bits<20; bits++)
555 {
556 io_length=*length;
557 packed=Ptngc_pack_array(coder,input,&io_length,
558 TNG_COMPRESS_ALGO_STOPBIT1,bits,natoms,0);
559 if (packed)
560 {
561 if ((new_parameter==-1) || (io_length<best_length))
562 {
563 new_parameter=bits;
564 best_length=io_length;
565 }
566 free(packed);
567 }
568 }
569 if (new_parameter==-1)
570 return 1;
571
572 *coding_parameter=new_parameter;
573 *length=best_length;
574 return 0;
575}
576
577static int determine_best_coding_triple(struct coder *coder,int *input, int *length,
578 int *coding_parameter, int natoms)
579{
580 int bits;
581 unsigned char *packed;
582 int best_length=0;
583 int new_parameter=-1;
584 int io_length;
585 for (bits=1; bits<20; bits++)
586 {
587 io_length=*length;
588 packed=Ptngc_pack_array(coder,input,&io_length,
589 TNG_COMPRESS_ALGO_TRIPLET2,bits,natoms,0);
590 if (packed)
591 {
592 if ((new_parameter==-1) || (io_length<best_length))
593 {
594 new_parameter=bits;
595 best_length=io_length;
596 }
597 free(packed);
598 }
599 }
600 if (new_parameter==-1)
601 return 1;
602
603 *coding_parameter=new_parameter;
604 *length=best_length;
605 return 0;
606}
607
608static void determine_best_pos_initial_coding(int *quant, int *quant_intra, int natoms, int speed,
609 fix_t prec_hi, fix_t prec_lo,
610 int *initial_coding, int *initial_coding_parameter)
611{
612 if (*initial_coding==-1)
613 {
614 /* Determine all parameters automatically */
615 int best_coding;
616 int best_coding_parameter;
617 int best_code_size;
618 int current_coding;
619 int current_coding_parameter;
620 int current_code_size;
621 struct coder *coder;
622 /* Start with XTC2, it should always work. */
623 current_coding=TNG_COMPRESS_ALGO_POS_XTC25;
624 current_coding_parameter=0;
625 compress_quantized_pos(quant,NULL((void*)0),quant_intra,natoms,1,speed,
626 current_coding,current_coding_parameter,
627 0,0,prec_hi,prec_lo,&current_code_size,NULL((void*)0));
628 best_coding=current_coding;
629 best_coding_parameter=current_coding_parameter;
630 best_code_size=current_code_size;
631
632 /* Determine best parameter for triplet intra. */
633 current_coding=TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3;
634 coder=Ptngc_coder_init();
635 current_code_size=natoms*3;
636 current_coding_parameter=0;
637 if (!determine_best_coding_triple(coder,quant_intra,&current_code_size,&current_coding_parameter,natoms))
638 {
639 if (current_code_size<best_code_size)
640 {
641 best_coding=current_coding;
642 best_coding_parameter=current_coding_parameter;
643 best_code_size=current_code_size;
644 }
645 }
646 Ptngc_coder_deinit(coder);
647
648 /* Determine best parameter for triplet one-to-one. */
649 current_coding=TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7;
650 coder=Ptngc_coder_init();
651 current_code_size=natoms*3;
652 current_coding_parameter=0;
653 if (!determine_best_coding_triple(coder,quant,&current_code_size,&current_coding_parameter,natoms))
654 {
655 if (current_code_size<best_code_size)
656 {
657 best_coding=current_coding;
658 best_coding_parameter=current_coding_parameter;
659 best_code_size=current_code_size;
660 }
661 }
662 Ptngc_coder_deinit(coder);
663
664 if (speed>=2)
665 {
666 current_coding=TNG_COMPRESS_ALGO_POS_XTC310;
667 current_coding_parameter=0;
668 compress_quantized_pos(quant,NULL((void*)0),quant_intra,natoms,1,speed,
669 current_coding,current_coding_parameter,
670 0,0,prec_hi,prec_lo,&current_code_size,NULL((void*)0));
671 if (current_code_size<best_code_size)
672 {
673 best_coding=current_coding;
674 best_coding_parameter=current_coding_parameter;
675 best_code_size=current_code_size;
676 }
677 }
678 /* Test BWLZH intra */
679 if (speed>=6)
680 {
681 current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9;
682 current_coding_parameter=0;
683 compress_quantized_pos(quant,NULL((void*)0),quant_intra,natoms,1,speed,
684 current_coding,current_coding_parameter,
685 0,0,prec_hi,prec_lo,&current_code_size,NULL((void*)0));
686 if (current_code_size<best_code_size)
687 {
688 best_coding=current_coding;
689 best_coding_parameter=current_coding_parameter;
690 }
691 }
692 *initial_coding=best_coding;
693 *initial_coding_parameter=best_coding_parameter;
694 }
695 else
696 {
697 if (*initial_coding_parameter==-1)
698 {
699 if ((*initial_coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
700 (*initial_coding==TNG_COMPRESS_ALGO_POS_XTC310) ||
701 (*initial_coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
702 *initial_coding_parameter=0;
703 else if (*initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3)
704 {
705 struct coder *coder=Ptngc_coder_init();
706 int current_code_size=natoms*3;
707 determine_best_coding_triple(coder,quant_intra,&current_code_size,initial_coding_parameter,natoms);
708 Ptngc_coder_deinit(coder);
709 }
710 else if (*initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7)
711 {
712 struct coder *coder=Ptngc_coder_init();
713 int current_code_size=natoms*3;
714 determine_best_coding_triple(coder,quant,&current_code_size,initial_coding_parameter,natoms);
715 Ptngc_coder_deinit(coder);
716 }
717 }
718 }
719}
720
721static void determine_best_pos_coding(int *quant, int *quant_inter, int *quant_intra, int natoms, int nframes, int speed,
722 fix_t prec_hi, fix_t prec_lo,
723 int *coding, int *coding_parameter)
724{
725 if (*coding==-1)
726 {
727 /* Determine all parameters automatically */
728 int best_coding;
729 int best_coding_parameter;
730 int best_code_size;
731 int current_coding;
732 int current_coding_parameter;
733 int current_code_size;
734 int initial_code_size;
735 struct coder *coder;
736 /* Always use XTC2 for the initial coding. */
737 compress_quantized_pos(quant,quant_inter,quant_intra,natoms,1,speed,
738 TNG_COMPRESS_ALGO_POS_XTC25,0,
739 0,0,
740 prec_hi,prec_lo,&initial_code_size,NULL((void*)0));
741 /* Start with XTC2, it should always work. */
742 current_coding=TNG_COMPRESS_ALGO_POS_XTC25;
743 current_coding_parameter=0;
744 compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed,
745 TNG_COMPRESS_ALGO_POS_XTC25,0,
746 current_coding,current_coding_parameter,
747 prec_hi,prec_lo,&current_code_size,NULL((void*)0));
748 best_coding=current_coding;
749 best_coding_parameter=current_coding_parameter;
750 best_code_size=current_code_size-initial_code_size; /* Correct for the use of XTC2 for the first frame. */
751
752 /* Determine best parameter for stopbit interframe coding. */
753 current_coding=TNG_COMPRESS_ALGO_POS_STOPBIT_INTER1;
754 coder=Ptngc_coder_init();
755 current_code_size=natoms*3*(nframes-1);
756 current_coding_parameter=0;
757 if (!determine_best_coding_stop_bits(coder,quant_inter+natoms*3,&current_code_size,
758 &current_coding_parameter,natoms))
759 {
760 if (current_code_size<best_code_size)
761 {
762 best_coding=current_coding;
763 best_coding_parameter=current_coding_parameter;
764 best_code_size=current_code_size;
765 }
766 }
767 Ptngc_coder_deinit(coder);
768
769 /* Determine best parameter for triplet interframe coding. */
770 current_coding=TNG_COMPRESS_ALGO_POS_TRIPLET_INTER2;
771 coder=Ptngc_coder_init();
772 current_code_size=natoms*3*(nframes-1);
773 current_coding_parameter=0;
774 if (!determine_best_coding_triple(coder,quant_inter+natoms*3,&current_code_size,
775 &current_coding_parameter,natoms))
776 {
777 if (current_code_size<best_code_size)
778 {
779 best_coding=current_coding;
780 best_coding_parameter=current_coding_parameter;
781 best_code_size=current_code_size;
782 }
783 }
784 Ptngc_coder_deinit(coder);
785
786 /* Determine best parameter for triplet intraframe coding. */
787 current_coding=TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3;
788 coder=Ptngc_coder_init();
789 current_code_size=natoms*3*(nframes-1);
790 current_coding_parameter=0;
791 if (!determine_best_coding_triple(coder,quant_intra+natoms*3,&current_code_size,
792 &current_coding_parameter,natoms))
793 {
794 if (current_code_size<best_code_size)
795 {
796 best_coding=current_coding;
797 best_coding_parameter=current_coding_parameter;
798 best_code_size=current_code_size;
799 }
800 }
801 Ptngc_coder_deinit(coder);
802
803 /* Determine best parameter for triplet one-to-one coding. */
804 current_coding=TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7;
805 coder=Ptngc_coder_init();
806 current_code_size=natoms*3*(nframes-1);
807 current_coding_parameter=0;
808 if (!determine_best_coding_triple(coder,quant+natoms*3,&current_code_size,
809 &current_coding_parameter,natoms))
810 {
811 if (current_code_size<best_code_size)
812 {
813 best_coding=current_coding;
814 best_coding_parameter=current_coding_parameter;
815 best_code_size=current_code_size;
816 }
817 }
818 Ptngc_coder_deinit(coder);
819
820 /* Test BWLZH inter */
821 if (speed>=4)
822 {
823 current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTER8;
824 current_coding_parameter=0;
825 compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed,
826 TNG_COMPRESS_ALGO_POS_XTC25,0,
827 current_coding,current_coding_parameter,
828 prec_hi,prec_lo,&current_code_size,NULL((void*)0));
829 current_code_size-=initial_code_size; /* Correct for the use of XTC2 for the first frame. */
830 if (current_code_size<best_code_size)
831 {
832 best_coding=current_coding;
833 best_coding_parameter=current_coding_parameter;
834 best_code_size=current_code_size;
835 }
836 }
837
838 /* Test BWLZH intra */
839 if (speed>=6)
840 {
841 current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9;
842 current_coding_parameter=0;
843 compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed,
844 TNG_COMPRESS_ALGO_POS_XTC25,0,
845 current_coding,current_coding_parameter,
846 prec_hi,prec_lo,&current_code_size,NULL((void*)0));
847 current_code_size-=initial_code_size; /* Correct for the use of XTC2 for the first frame. */
848 if (current_code_size<best_code_size)
849 {
850 best_coding=current_coding;
851 best_coding_parameter=current_coding_parameter;
852 }
853 }
854 *coding=best_coding;
855 *coding_parameter=best_coding_parameter;
856 }
857 else if (*coding_parameter==-1)
858 {
859 if ((*coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
860 (*coding==TNG_COMPRESS_ALGO_POS_XTC310) ||
861 (*coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTER8) ||
862 (*coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
863 *coding_parameter=0;
864 else if (*coding==TNG_COMPRESS_ALGO_POS_STOPBIT_INTER1)
865 {
866 struct coder *coder=Ptngc_coder_init();
867 int current_code_size=natoms*3*(nframes-1);
868 determine_best_coding_stop_bits(coder,quant_inter+natoms*3,&current_code_size,
869 coding_parameter,natoms);
870 Ptngc_coder_deinit(coder);
871 }
872 else if (*coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTER2)
873 {
874 struct coder *coder=Ptngc_coder_init();
875 int current_code_size=natoms*3*(nframes-1);
876 determine_best_coding_triple(coder,quant_inter+natoms*3,&current_code_size,
877 coding_parameter,natoms);
878 Ptngc_coder_deinit(coder);
879 }
880 else if (*coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3)
881 {
882 struct coder *coder=Ptngc_coder_init();
883 int current_code_size=natoms*3*(nframes-1);
884 determine_best_coding_triple(coder,quant_intra+natoms*3,&current_code_size,
885 coding_parameter,natoms);
886 Ptngc_coder_deinit(coder);
887 }
888 else if (*coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7)
889 {
890 struct coder *coder=Ptngc_coder_init();
891 int current_code_size=natoms*3*(nframes-1);
892 determine_best_coding_triple(coder,quant+natoms*3,&current_code_size,
893 coding_parameter,natoms);
894 Ptngc_coder_deinit(coder);
895 }
896 }
897}
898
899static void determine_best_vel_initial_coding(int *quant, int natoms, int speed,
900 fix_t prec_hi, fix_t prec_lo,
901 int *initial_coding, int *initial_coding_parameter)
902{
903 if (*initial_coding==-1)
904 {
905 /* Determine all parameters automatically */
906 int best_coding=-1;
907 int best_coding_parameter=-1;
908 int best_code_size=-1;
909 int current_coding;
910 int current_coding_parameter;
911 int current_code_size;
912 struct coder *coder;
913 /* Start to determine best parameter for stopbit one-to-one. */
914 current_coding=TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1;
915 current_code_size=natoms*3;
916 current_coding_parameter=0;
917 coder=Ptngc_coder_init();
918 if (!determine_best_coding_stop_bits(coder,quant,&current_code_size,
919 &current_coding_parameter,natoms))
920 {
921 best_coding=current_coding;
922 best_coding_parameter=current_coding_parameter;
923 best_code_size=current_code_size;
924 }
925 Ptngc_coder_deinit(coder);
926
927 /* Determine best parameter for triplet one-to-one. */
928 current_coding=TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3;
929 coder=Ptngc_coder_init();
930 current_code_size=natoms*3;
931 current_coding_parameter=0;
932 if (!determine_best_coding_triple(coder,quant,&current_code_size,&current_coding_parameter,natoms))
933 {
934 if ((best_coding==-1) || (current_code_size<best_code_size))
935 {
936 best_coding=current_coding;
937 best_coding_parameter=current_coding_parameter;
938 best_code_size=current_code_size;
939 }
940 }
941 Ptngc_coder_deinit(coder);
942
943 /* Test BWLZH one-to-one */
944 if (speed>=4)
945 {
946 current_coding=TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9;
947 current_coding_parameter=0;
948 compress_quantized_vel(quant,NULL((void*)0),natoms,1,speed,
949 current_coding,current_coding_parameter,
950 0,0,prec_hi,prec_lo,&current_code_size,NULL((void*)0));
951 if ((best_coding==-1) || (current_code_size<best_code_size))
952 {
953 best_coding=current_coding;
954 best_coding_parameter=current_coding_parameter;
955 }
956 }
957 *initial_coding=best_coding;
958 *initial_coding_parameter=best_coding_parameter;
959 }
960 else if (*initial_coding_parameter==-1)
961 {
962 if (*initial_coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9)
963 *initial_coding_parameter=0;
964 else if (*initial_coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1)
965 {
966 struct coder *coder=Ptngc_coder_init();
967 int current_code_size=natoms*3;
968 determine_best_coding_stop_bits(coder,quant,&current_code_size,
969 initial_coding_parameter,natoms);
970 Ptngc_coder_deinit(coder);
971 }
972 else if (*initial_coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3)
973 {
974 struct coder *coder=Ptngc_coder_init();
975 int current_code_size=natoms*3;
976 determine_best_coding_triple(coder,quant,&current_code_size,initial_coding_parameter,natoms);
977 Ptngc_coder_deinit(coder);
978 }
979 }
980}
981
982static void determine_best_vel_coding(int *quant, int *quant_inter, int natoms, int nframes, int speed,
983 fix_t prec_hi, fix_t prec_lo,
984 int *coding, int *coding_parameter)
985{
986 if (*coding==-1)
987 {
988 /* Determine all parameters automatically */
989 int best_coding;
990 int best_coding_parameter;
991 int best_code_size;
992 int current_coding;
993 int current_coding_parameter;
994 int current_code_size;
995 int initial_code_size;
996 int initial_numbits=5;
997 struct coder *coder;
998 /* Use stopbits one-to-one coding for the initial coding. */
999 compress_quantized_vel(quant,NULL((void*)0),natoms,1,speed,
1000 TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1,initial_numbits,
1001 0,0,prec_hi,prec_lo,&initial_code_size,NULL((void*)0));
1002
1003 /* Test stopbit one-to-one */
1004 current_coding=TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1;
1005 current_code_size=natoms*3*(nframes-1);
1006 current_coding_parameter=0;
1007 coder=Ptngc_coder_init();
1008 determine_best_coding_stop_bits(coder,quant+natoms*3,&current_code_size,
1009 &current_coding_parameter,natoms);
1010 Ptngc_coder_deinit(coder);
1011 best_coding=current_coding;
1012 best_code_size=current_code_size;
1013 best_coding_parameter=current_coding_parameter;
1014
1015 /* Test triplet interframe */
1016 current_coding=TNG_COMPRESS_ALGO_VEL_TRIPLET_INTER2;
1017 current_code_size=natoms*3*(nframes-1);
1018 current_coding_parameter=0;
1019 coder=Ptngc_coder_init();
1020 if (!determine_best_coding_triple(coder,quant_inter+natoms*3,&current_code_size,
1021 &current_coding_parameter,natoms))
1022 {
1023 if (current_code_size<best_code_size)
1024 {
1025 best_coding=current_coding;
1026 best_code_size=current_code_size;
1027 best_coding_parameter=current_coding_parameter;
1028 }
1029 }
1030 Ptngc_coder_deinit(coder);
1031
1032 /* Test triplet one-to-one */
1033 current_coding=TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3;
1034 current_code_size=natoms*3*(nframes-1);
1035 current_coding_parameter=0;
1036 coder=Ptngc_coder_init();
1037 if (!determine_best_coding_triple(coder,quant+natoms*3,&current_code_size,
1038 &current_coding_parameter,natoms))
1039 {
1040 if (current_code_size<best_code_size)
1041 {
1042 best_coding=current_coding;
1043 best_code_size=current_code_size;
1044 best_coding_parameter=current_coding_parameter;
1045 }
1046 }
1047 Ptngc_coder_deinit(coder);
1048
1049 /* Test stopbit interframe */
1050 current_coding=TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER6;
1051 current_code_size=natoms*3*(nframes-1);
1052 current_coding_parameter=0;
1053 coder=Ptngc_coder_init();
1054 if (!determine_best_coding_stop_bits(coder,quant_inter+natoms*3,&current_code_size,
1055 &current_coding_parameter,natoms))
1056 {
1057 if (current_code_size<best_code_size)
1058 {
1059 best_coding=current_coding;
1060 best_code_size=current_code_size;
1061 best_coding_parameter=current_coding_parameter;
1062 }
1063 }
1064 Ptngc_coder_deinit(coder);
1065
1066 if (speed>=4)
1067 {
1068 /* Test BWLZH inter */
1069 current_coding=TNG_COMPRESS_ALGO_VEL_BWLZH_INTER8;
1070 current_coding_parameter=0;
1071 compress_quantized_vel(quant,quant_inter,natoms,nframes,speed,
1072 TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1,initial_numbits,
1073 current_coding,current_coding_parameter,
1074 prec_hi,prec_lo,&current_code_size,NULL((void*)0));
1075 current_code_size-=initial_code_size; /* Correct for the initial frame */
1076 if (current_code_size<best_code_size)
1077 {
1078 best_coding=current_coding;
1079 best_code_size=current_code_size;
1080 best_coding_parameter=current_coding_parameter;
1081 }
1082
1083 /* Test BWLZH one-to-one */
1084 current_coding=TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9;
1085 current_coding_parameter=0;
1086 compress_quantized_vel(quant,quant_inter,natoms,nframes,speed,
1087 TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1,initial_numbits,
1088 current_coding,current_coding_parameter,
1089 prec_hi,prec_lo,&current_code_size,NULL((void*)0));
1090 current_code_size-=initial_code_size; /* Correct for the initial frame */
1091 if (current_code_size<best_code_size)
1092 {
1093 best_coding=current_coding;
1094 best_coding_parameter=current_coding_parameter;
1095 }
1096 }
1097 *coding=best_coding;
1098 *coding_parameter=best_coding_parameter;
1099 }
1100 else if (*coding_parameter==-1)
1101 {
1102 if ((*coding==TNG_COMPRESS_ALGO_VEL_BWLZH_INTER8) ||
1103 (*coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9))
1104 *coding_parameter=0;
1105 else if (*coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1)
1106 {
1107 struct coder *coder=Ptngc_coder_init();
1108 int current_code_size=natoms*3*(nframes-1);
1109 determine_best_coding_stop_bits(coder,quant+natoms*3,&current_code_size,
1110 coding_parameter,natoms);
1111 Ptngc_coder_deinit(coder);
1112 }
1113 else if (*coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_INTER2)
1114 {
1115 struct coder *coder=Ptngc_coder_init();
1116 int current_code_size=natoms*3*(nframes-1);
1117 determine_best_coding_triple(coder,quant_inter+natoms*3,&current_code_size,
1118 coding_parameter,natoms);
1119 Ptngc_coder_deinit(coder);
1120 }
1121 else if (*coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3)
1122 {
1123 struct coder *coder=Ptngc_coder_init();
1124 int current_code_size=natoms*3*(nframes-1);
1125 determine_best_coding_triple(coder,quant+natoms*3,&current_code_size,
1126 coding_parameter,natoms);
1127 Ptngc_coder_deinit(coder);
1128 }
1129 else if (*coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER6)
1130 {
1131 struct coder *coder=Ptngc_coder_init();
1132 int current_code_size=natoms*3*(nframes-1);
1133 determine_best_coding_stop_bits(coder,quant_inter+natoms*3,&current_code_size,
1134 coding_parameter,natoms);
1135 Ptngc_coder_deinit(coder);
1136 }
1137 }
1138}
1139
1140char DECLSPECDLLEXPORT *tng_compress_pos_int(int *pos, int natoms, int nframes,
1141 unsigned long prec_hi, unsigned long prec_lo,
1142 int speed,int *algo,
1143 int *nitems)
1144{
1145 char *data=malloc(natoms*nframes*14+11*4); /* 12 bytes are required to store 4 32 bit integers
1146 This is 17% extra. The final 11*4 is to store information
1147 needed for decompression. */
1148 int *quant=pos; /* Already quantized positions. */
1149 int *quant_intra=malloc(natoms*nframes*3*sizeof *quant_intra);
1150 int *quant_inter=malloc(natoms*nframes*3*sizeof *quant_inter);
1151
1152 int initial_coding, initial_coding_parameter;
1153 int coding, coding_parameter;
1154 if (speed==0)
1155 speed=SPEED_DEFAULT2; /* Set the default speed. */
1156 /* Boundaries of speed. */
1157 if (speed<1)
1158 speed=1;
1159 if (speed>6)
1160 speed=6;
1161 initial_coding=algo[0];
1162 initial_coding_parameter=algo[1];
1163 coding=algo[2];
1164 coding_parameter=algo[3];
1165
1166 quant_inter_differences(quant,natoms,nframes,quant_inter);
1167 quant_intra_differences(quant,natoms,nframes,quant_intra);
1168
1169 /* If any of the above codings / coding parameters are == -1, the optimal parameters must be found */
1170 if (initial_coding==-1)
1171 {
1172 initial_coding_parameter=-1;
1173 determine_best_pos_initial_coding(quant,quant_intra,natoms,speed,prec_hi,prec_lo,
1174 &initial_coding,&initial_coding_parameter);
1175 }
1176 else if (initial_coding_parameter==-1)
1177 {
1178 determine_best_pos_initial_coding(quant,quant_intra,natoms,speed,prec_hi,prec_lo,
1179 &initial_coding,&initial_coding_parameter);
1180 }
1181
1182 if (nframes==1)
1183 {
1184 coding=0;
1185 coding_parameter=0;
1186 }
1187
1188 if (nframes>1)
1189 {
1190 if (coding==-1)
1191 {
1192 coding_parameter=-1;
1193 determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo,
1194 &coding,&coding_parameter);
1195 }
1196 else if (coding_parameter==-1)
1197 {
1198 determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo,
1199 &coding,&coding_parameter);
1200 }
1201 }
1202
1203 compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed,
1204 initial_coding,initial_coding_parameter,
1205 coding,coding_parameter,
1206 prec_hi,prec_lo,nitems,data);
1207 free(quant_inter);
1208 free(quant_intra);
1209 if (algo[0]==-1)
1210 algo[0]=initial_coding;
1211 if (algo[1]==-1)
1212 algo[1]=initial_coding_parameter;
1213 if (algo[2]==-1)
1214 algo[2]=coding;
1215 if (algo[3]==-1)
1216 algo[3]=coding_parameter;
1217 return data;
1218}
1219
1220char DECLSPECDLLEXPORT *tng_compress_pos(double *pos, int natoms, int nframes,
1221 double desired_precision,
1222 int speed,int *algo,
1223 int *nitems)
1224{
1225 int *quant=malloc(natoms*nframes*3*sizeof *quant);
1226 char *data;
1227 fix_t prec_hi, prec_lo;
1228 Ptngc_d_to_i32x2(desired_precision,&prec_hi,&prec_lo);
1229
1230 if (quantize(pos,natoms,nframes,PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),quant))
1231 data=NULL((void*)0); /* Error occured. Too large input values. */
1232 else
1233 data=tng_compress_pos_int(quant,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1234 free(quant);
1235 return data;
1236}
1237
1238char DECLSPECDLLEXPORT *tng_compress_pos_float(float *pos, int natoms, int nframes,
1239 float desired_precision,
1240 int speed,int *algo,
1241 int *nitems)
1242{
1243 int *quant=malloc(natoms*nframes*3*sizeof *quant);
1244 char *data;
1245 fix_t prec_hi, prec_lo;
1246 Ptngc_d_to_i32x2((double)desired_precision,&prec_hi,&prec_lo);
1247
1248 if (quantize_float(pos,natoms,nframes,(float)PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),quant))
1249 data=NULL((void*)0); /* Error occured. Too large input values. */
1250 else
1251 data=tng_compress_pos_int(quant,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1252 free(quant);
1253 return data;
1254}
1255
1256char DECLSPECDLLEXPORT *tng_compress_pos_find_algo(double *pos, int natoms, int nframes,
1257 double desired_precision,
1258 int speed,
1259 int *algo,
1260 int *nitems)
1261{
1262 algo[0]=-1;
1263 algo[1]=-1;
1264 algo[2]=-1;
1265 algo[3]=-1;
1266 return tng_compress_pos(pos,natoms,nframes,desired_precision,speed,algo,nitems);
1267}
1268
1269char DECLSPECDLLEXPORT *tng_compress_pos_float_find_algo(float *pos, int natoms, int nframes,
1270 float desired_precision,
1271 int speed,
1272 int *algo,
1273 int *nitems)
1274{
1275 algo[0]=-1;
1276 algo[1]=-1;
1277 algo[2]=-1;
1278 algo[3]=-1;
1279 return tng_compress_pos_float(pos,natoms,nframes,desired_precision,speed,algo,nitems);
1280}
1281
1282char DECLSPECDLLEXPORT *tng_compress_pos_int_find_algo(int *pos, int natoms, int nframes,
1283 unsigned long prec_hi, unsigned long prec_lo,
1284 int speed,int *algo,
1285 int *nitems)
1286{
1287 algo[0]=-1;
1288 algo[1]=-1;
1289 algo[2]=-1;
1290 algo[3]=-1;
1291 return tng_compress_pos_int(pos,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1292}
1293
1294
1295
1296int DECLSPECDLLEXPORT tng_compress_nalgo(void)
1297{
1298 return 4; /* There are currently four parameters required:
1299
1300 1) The compression algorithm for the first frame (initial_coding).
1301 2) One parameter to the algorithm for the first frame (the initial coding parameter).
1302 3) The compression algorithm for the remaining frames (coding).
1303 4) One parameter to the algorithm for the remaining frames (the coding parameter). */
1304}
1305
1306char DECLSPECDLLEXPORT *tng_compress_vel_int(int *vel, int natoms, int nframes,
1307 unsigned long prec_hi, unsigned long prec_lo,
1308 int speed, int *algo,
1309 int *nitems)
1310{
1311 char *data=malloc(natoms*nframes*14+11*4); /* 12 bytes are required to store 4 32 bit integers
1312 This is 17% extra. The final 11*4 is to store information
1313 needed for decompression. */
1314 int *quant=vel;
1315 int *quant_inter=malloc(natoms*nframes*3*sizeof *quant_inter);
1316
1317 int initial_coding, initial_coding_parameter;
1318 int coding, coding_parameter;
1319 if (speed==0)
1320 speed=SPEED_DEFAULT2; /* Set the default speed. */
1321 /* Boundaries of speed. */
1322 if (speed<1)
1323 speed=1;
1324 if (speed>6)
1325 speed=6;
1326 initial_coding=algo[0];
1327 initial_coding_parameter=algo[1];
1328 coding=algo[2];
1329 coding_parameter=algo[3];
1330
1331 quant_inter_differences(quant,natoms,nframes,quant_inter);
1332
1333 /* If any of the above codings / coding parameters are == -1, the optimal parameters must be found */
1334 if (initial_coding==-1)
1335 {
1336 initial_coding_parameter=-1;
1337 determine_best_vel_initial_coding(quant,natoms,speed,prec_hi,prec_lo,
1338 &initial_coding,&initial_coding_parameter);
1339 }
1340 else if (initial_coding_parameter==-1)
1341 {
1342 determine_best_vel_initial_coding(quant,natoms,speed,prec_hi,prec_lo,
1343 &initial_coding,&initial_coding_parameter);
1344 }
1345
1346 if (nframes==1)
1347 {
1348 coding=0;
1349 coding_parameter=0;
1350 }
1351
1352 if (nframes>1)
1353 {
1354 if (coding==-1)
1355 {
1356 coding_parameter=-1;
1357 determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo,
1358 &coding,&coding_parameter);
1359 }
1360 else if (coding_parameter==-1)
1361 {
1362 determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo,
1363 &coding,&coding_parameter);
1364 }
1365 }
1366
1367 compress_quantized_vel(quant,quant_inter,natoms,nframes,speed,
1368 initial_coding,initial_coding_parameter,
1369 coding,coding_parameter,
1370 prec_hi,prec_lo,nitems,data);
1371 free(quant_inter);
1372 if (algo[0]==-1)
1373 algo[0]=initial_coding;
1374 if (algo[1]==-1)
1375 algo[1]=initial_coding_parameter;
1376 if (algo[2]==-1)
1377 algo[2]=coding;
1378 if (algo[3]==-1)
1379 algo[3]=coding_parameter;
1380 return data;
1381}
1382
1383char DECLSPECDLLEXPORT *tng_compress_vel(double *vel, int natoms, int nframes,
1384 double desired_precision,
1385 int speed, int *algo,
1386 int *nitems)
1387{
1388 int *quant=malloc(natoms*nframes*3*sizeof *quant);
1389 char *data;
1390 fix_t prec_hi, prec_lo;
1391 Ptngc_d_to_i32x2(desired_precision,&prec_hi,&prec_lo);
1392 if (quantize(vel,natoms,nframes,PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),quant))
1393 data=NULL((void*)0); /* Error occured. Too large input values. */
1394 else
1395 data=tng_compress_vel_int(quant,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1396 free(quant);
1397 return data;
1398}
1399
1400char DECLSPECDLLEXPORT *tng_compress_vel_float(float *vel, int natoms, int nframes,
1401 float desired_precision,
1402 int speed, int *algo,
1403 int *nitems)
1404{
1405 int *quant=malloc(natoms*nframes*3*sizeof *quant);
1406 char *data;
1407 fix_t prec_hi, prec_lo;
1408 Ptngc_d_to_i32x2((double)desired_precision,&prec_hi,&prec_lo);
1409 if (quantize_float(vel,natoms,nframes,(float)PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),quant))
1410 data=NULL((void*)0); /* Error occured. Too large input values. */
1411 else
1412 data=tng_compress_vel_int(quant,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1413 free(quant);
1414 return data;
1415}
1416
1417char DECLSPECDLLEXPORT *tng_compress_vel_find_algo(double *vel, int natoms, int nframes,
1418 double desired_precision,
1419 int speed,
1420 int *algo,
1421 int *nitems)
1422{
1423 algo[0]=-1;
1424 algo[1]=-1;
1425 algo[2]=-1;
1426 algo[3]=-1;
1427 return tng_compress_vel(vel,natoms,nframes,desired_precision,speed,algo,nitems);
1428}
1429
1430char DECLSPECDLLEXPORT *tng_compress_vel_float_find_algo(float *vel, int natoms, int nframes,
1431 float desired_precision,
1432 int speed,
1433 int *algo,
1434 int *nitems)
1435{
1436 algo[0]=-1;
1437 algo[1]=-1;
1438 algo[2]=-1;
1439 algo[3]=-1;
1440 return tng_compress_vel_float(vel,natoms,nframes,desired_precision,speed,algo,nitems);
1441}
1442
1443char DECLSPECDLLEXPORT *tng_compress_vel_int_find_algo(int *vel, int natoms, int nframes,
1444 unsigned long prec_hi, unsigned long prec_lo,
1445 int speed,
1446 int *algo,
1447 int *nitems)
1448{
1449 algo[0]=-1;
1450 algo[1]=-1;
1451 algo[2]=-1;
1452 algo[3]=-1;
1453 return tng_compress_vel_int(vel,natoms,nframes,prec_hi,prec_lo,speed,algo,nitems);
1454}
1455
1456int DECLSPECDLLEXPORT tng_compress_inquire(char *data,int *vel, int *natoms,
1457 int *nframes, double *precision,
1458 int *algo)
1459{
1460 int bufloc=0;
1461 fix_t prec_hi, prec_lo;
1462 int initial_coding, initial_coding_parameter;
1463 int coding, coding_parameter;
1464 int magic_int;
1465 magic_int=(int)readbufferfix((unsigned char *)data+bufloc,4);
1466 bufloc+=4;
1467 if (magic_int==MAGIC_INT_POS0x50474E54)
1468 *vel=0;
1469 else if (magic_int==MAGIC_INT_VEL0x56474E54)
1470 *vel=1;
1471 else
1472 return 1;
1473 /* Number of atoms. */
1474 *natoms=(int)readbufferfix((unsigned char *)data+bufloc,4);
1475 bufloc+=4;
1476 /* Number of frames. */
1477 *nframes=(int)readbufferfix((unsigned char *)data+bufloc,4);
1478 bufloc+=4;
1479 /* Initial coding. */
1480 initial_coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1481 bufloc+=4;
1482 /* Initial coding parameter. */
1483 initial_coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1484 bufloc+=4;
1485 /* Coding. */
1486 coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1487 bufloc+=4;
1488 /* Coding parameter. */
1489 coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1490 bufloc+=4;
1491 /* Precision. */
1492 prec_lo=readbufferfix((unsigned char *)data+bufloc,4);
1493 bufloc+=4;
1494 prec_hi=readbufferfix((unsigned char *)data+bufloc,4);
1495 *precision=PRECISION(prec_hi, prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo));
1496 algo[0]=initial_coding;
1497 algo[1]=initial_coding_parameter;
1498 algo[2]=coding;
1499 algo[3]=coding_parameter;
1500 return 0;
1501}
1502
1503static int tng_compress_uncompress_pos_gen(char *data,double *posd,float *posf,int *posi,unsigned long *prec_hi, unsigned long *prec_lo)
1504{
1505 int bufloc=0;
1506 int length;
1507 int natoms, nframes;
1508 int initial_coding, initial_coding_parameter;
1509 int coding, coding_parameter;
1510 int *quant=NULL((void*)0);
1511 struct coder *coder=NULL((void*)0);
1512 int rval=0;
1513 int magic_int;
1514 /* Magic integer for positions */
1515 magic_int=(int)readbufferfix((unsigned char *)data+bufloc,4);
1516 bufloc+=4;
1517 if (magic_int!=MAGIC_INT_POS0x50474E54)
1518 {
1519 rval=1;
1520 goto error;
1521 }
1522 /* Number of atoms. */
1523 natoms=(int)readbufferfix((unsigned char *)data+bufloc,4);
1524 bufloc+=4;
1525 /* Number of frames. */
1526 nframes=(int)readbufferfix((unsigned char *)data+bufloc,4);
1527 bufloc+=4;
1528 /* Initial coding. */
1529 initial_coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1530 bufloc+=4;
1531 /* Initial coding parameter. */
1532 initial_coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1533 bufloc+=4;
1534 /* Coding. */
1535 coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1536 bufloc+=4;
1537 /* Coding parameter. */
1538 coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1539 bufloc+=4;
1540 /* Precision. */
1541 *prec_lo=readbufferfix((unsigned char *)data+bufloc,4);
1542 bufloc+=4;
1543 *prec_hi=readbufferfix((unsigned char *)data+bufloc,4);
1544 bufloc+=4;
1545 /* Allocate the memory for the quantized positions */
1546 quant=malloc(natoms*nframes*3*sizeof *quant);
1547 /* The data block length. */
1548 length=(int)readbufferfix((unsigned char *)data+bufloc,4);
1549 bufloc+=4;
1550 /* The initial frame */
1551 coder=Ptngc_coder_init();
1552 rval=Ptngc_unpack_array(coder,(unsigned char*)data+bufloc,quant,natoms*3,
1553 initial_coding,initial_coding_parameter,natoms);
1554 Ptngc_coder_deinit(coder);
1555 if (rval)
1556 goto error;
1557 /* Skip past the actual data block. */
1558 bufloc+=length;
1559 /* Obtain the actual positions for the initial block. */
1560 if ((initial_coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
1561 (initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7) ||
1562 (initial_coding==TNG_COMPRESS_ALGO_POS_XTC310))
1563 {
1564 if (posd)
1565 unquantize(posd,natoms,1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1566 else if (posf)
1567 unquantize_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1568 else if (posi)
1569 memcpy(posi,quant,natoms*3*sizeof *posi);
1570 }
1571 else if ((initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3) ||
1572 (initial_coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
1573 {
1574 if (posd)
1575 unquantize_intra_differences(posd,natoms,1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1576 else if (posf)
1577 unquantize_intra_differences_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1578 else if (posi)
1579 unquantize_intra_differences_int(posi,natoms,1,quant);
1580 unquant_intra_differences_first_frame(quant,natoms);
1581 }
1582 /* The remaining frames. */
1583 if (nframes>1)
1584 {
1585 bufloc+=4;
1586 coder=Ptngc_coder_init();
1587 rval=Ptngc_unpack_array(coder,(unsigned char *)data+bufloc,quant+natoms*3,(nframes-1)*natoms*3,
1588 coding,coding_parameter,natoms);
1589 Ptngc_coder_deinit(coder);
1590 if (rval)
1591 goto error;
1592 if ((coding==TNG_COMPRESS_ALGO_POS_STOPBIT_INTER1) ||
1593 (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTER2) ||
1594 (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTER8))
1595 {
1596 /* This requires that the first frame is already in one-to-one format, even if intra-frame
1597 compression was done there. Therefore the unquant_intra_differences_first_frame should be called
1598 before to convert it correctly. */
1599 if (posd)
1600 unquantize_inter_differences(posd,natoms,nframes,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1601 else if (posf)
1602 unquantize_inter_differences_float(posf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1603 else if (posi)
1604 unquantize_inter_differences_int(posi,natoms,nframes,quant);
1605 }
1606 else if ((coding==TNG_COMPRESS_ALGO_POS_XTC25) ||
1607 (coding==TNG_COMPRESS_ALGO_POS_XTC310) ||
1608 (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE7))
1609 {
1610 if (posd)
1611 unquantize(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1612 else if (posf)
1613 unquantize_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1614 else if (posi)
1615 memcpy(posi+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *posi);
1616 }
1617 else if ((coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA3) ||
1618 (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA9))
1619 {
1620 if (posd)
1621 unquantize_intra_differences(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1622 else if (posf)
1623 unquantize_intra_differences_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1624 else if (posi)
1625 unquantize_intra_differences_int(posi+natoms*3,natoms,nframes-1,quant+natoms*3);
1626 }
1627 }
1628 error:
1629 free(quant);
1630 return rval;
1631}
1632
1633static int tng_compress_uncompress_pos(char *data,double *pos)
1634{
1635 unsigned long prec_hi, prec_lo;
1636 return tng_compress_uncompress_pos_gen(data,pos,NULL((void*)0),NULL((void*)0),&prec_hi,&prec_lo);
1637}
1638
1639static int tng_compress_uncompress_pos_float(char *data,float *pos)
1640{
1641 unsigned long prec_hi, prec_lo;
1642 return tng_compress_uncompress_pos_gen(data,NULL((void*)0),pos,NULL((void*)0),&prec_hi,&prec_lo);
1643}
1644
1645static int tng_compress_uncompress_pos_int(char *data,int *pos, unsigned long *prec_hi, unsigned long *prec_lo)
1646{
1647 return tng_compress_uncompress_pos_gen(data,NULL((void*)0),NULL((void*)0),pos,prec_hi,prec_lo);
1648}
1649
1650static int tng_compress_uncompress_vel_gen(char *data,double *veld,float *velf,int *veli,unsigned long *prec_hi, unsigned long *prec_lo)
1651{
1652 int bufloc=0;
1653 int length;
1654 int natoms, nframes;
1655 int initial_coding, initial_coding_parameter;
1656 int coding, coding_parameter;
1657 int *quant=NULL((void*)0);
1658 struct coder *coder=NULL((void*)0);
1659 int rval=0;
1660 int magic_int;
1661 /* Magic integer for velocities */
1662 magic_int=(int)readbufferfix((unsigned char *)data+bufloc,4);
1663 bufloc+=4;
1664 if (magic_int!=MAGIC_INT_VEL0x56474E54)
1665 {
1666 rval=1;
1667 goto error;
1668 }
1669 /* Number of atoms. */
1670 natoms=(int)readbufferfix((unsigned char *)data+bufloc,4);
1671 bufloc+=4;
1672 /* Number of frames. */
1673 nframes=(int)readbufferfix((unsigned char *)data+bufloc,4);
1674 bufloc+=4;
1675 /* Initial coding. */
1676 initial_coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1677 bufloc+=4;
1678 /* Initial coding parameter. */
1679 initial_coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1680 bufloc+=4;
1681 /* Coding. */
1682 coding=(int)readbufferfix((unsigned char *)data+bufloc,4);
1683 bufloc+=4;
1684 /* Coding parameter. */
1685 coding_parameter=(int)readbufferfix((unsigned char *)data+bufloc,4);
1686 bufloc+=4;
1687 /* Precision. */
1688 *prec_lo=readbufferfix((unsigned char *)data+bufloc,4);
1689 bufloc+=4;
1690 *prec_hi=readbufferfix((unsigned char *)data+bufloc,4);
1691 bufloc+=4;
1692 /* Allocate the memory for the quantized positions */
1693 quant=malloc(natoms*nframes*3*sizeof *quant);
1694 /* The data block length. */
1695 length=(int)readbufferfix((unsigned char *)data+bufloc,4);
1696 bufloc+=4;
1697 /* The initial frame */
1698 coder=Ptngc_coder_init();
1699 rval=Ptngc_unpack_array(coder,(unsigned char*)data+bufloc,quant,natoms*3,
1700 initial_coding,initial_coding_parameter,natoms);
1701 Ptngc_coder_deinit(coder);
1702 if (rval)
1703 goto error;
1704 /* Skip past the actual data block. */
1705 bufloc+=length;
1706 /* Obtain the actual positions for the initial block. */
1707 if ((initial_coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1) ||
1708 (initial_coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3) ||
1709 (initial_coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9))
1710 {
1711 if (veld)
1712 unquantize(veld,natoms,1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1713 else if (velf)
1714 unquantize_float(velf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1715 else if (veli)
1716 memcpy(veli,quant,natoms*3*sizeof *veli);
1717 }
1718 /* The remaining frames. */
1719 if (nframes>1)
1720 {
1721 bufloc+=4;
1722 coder=Ptngc_coder_init();
1723 rval=Ptngc_unpack_array(coder,(unsigned char *)data+bufloc,quant+natoms*3,(nframes-1)*natoms*3,
1724 coding,coding_parameter,natoms);
1725 Ptngc_coder_deinit(coder);
1726 if (rval)
1727 goto error;
1728 /* Inter-frame compression? */
1729 if ((coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_INTER2) ||
1730 (coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER6) ||
1731 (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_INTER8))
1732 {
1733 /* This requires that the first frame is already in one-to-one format. */
1734 if (veld)
1735 unquantize_inter_differences(veld,natoms,nframes,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1736 else if (velf)
1737 unquantize_inter_differences_float(velf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant);
1738 else if (veli)
1739 unquantize_inter_differences_int(veli,natoms,nframes,quant);
1740 }
1741 /* One-to-one compression? */
1742 else if ((coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE1) ||
1743 (coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE3) ||
1744 (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE9))
1745 {
1746 if (veld)
1747 unquantize(veld+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1748 else if (velf)
1749 unquantize_float(velf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo)(Ptngc_i32x2_to_d(*prec_hi,*prec_lo)),quant+natoms*3);
1750 else if (veli)
1751 memcpy(veli+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *veli);
1752 }
1753 }
1754 error:
1755 free(quant);
1756 return rval;
1757}
1758
1759static int tng_compress_uncompress_vel(char *data,double *vel)
1760{
1761 unsigned long prec_hi, prec_lo;
1762 return tng_compress_uncompress_vel_gen(data,vel,NULL((void*)0),NULL((void*)0),&prec_hi,&prec_lo);
1763}
1764
1765static int tng_compress_uncompress_vel_float(char *data,float *vel)
1766{
1767 unsigned long prec_hi, prec_lo;
1768 return tng_compress_uncompress_vel_gen(data,NULL((void*)0),vel,NULL((void*)0),&prec_hi,&prec_lo);
1769}
1770
1771static int tng_compress_uncompress_vel_int(char *data,int *vel, unsigned long *prec_hi, unsigned long *prec_lo)
1772{
1773 return tng_compress_uncompress_vel_gen(data,NULL((void*)0),NULL((void*)0),vel,prec_hi,prec_lo);
1774}
1775
1776/* Uncompresses any tng compress block, positions or velocities. It determines whether it is positions or velocities from the data buffer. The return value is 0 if ok, and 1 if not.
1777*/
1778int DECLSPECDLLEXPORT tng_compress_uncompress(char *data,double *posvel)
1779{
1780 int magic_int;
1781 magic_int=(int)readbufferfix((unsigned char *)data,4);
1782 if (magic_int==MAGIC_INT_POS0x50474E54)
1783 return tng_compress_uncompress_pos(data,posvel);
1784 else if (magic_int==MAGIC_INT_VEL0x56474E54)
1785 return tng_compress_uncompress_vel(data,posvel);
1786 else
1787 return 1;
1788}
1789
1790int DECLSPECDLLEXPORT tng_compress_uncompress_float(char *data,float *posvel)
1791{
1792 int magic_int;
1793 magic_int=(int)readbufferfix((unsigned char *)data,4);
1794 if (magic_int==MAGIC_INT_POS0x50474E54)
1795 return tng_compress_uncompress_pos_float(data,posvel);
1796 else if (magic_int==MAGIC_INT_VEL0x56474E54)
1797 return tng_compress_uncompress_vel_float(data,posvel);
1798 else
1799 return 1;
1800}
1801
1802int DECLSPECDLLEXPORT tng_compress_uncompress_int(char *data,int *posvel, unsigned long *prec_hi, unsigned long *prec_lo)
1803{
1804 int magic_int;
1805 magic_int=(int)readbufferfix((unsigned char *)data,4);
1806 if (magic_int==MAGIC_INT_POS0x50474E54)
1807 return tng_compress_uncompress_pos_int(data,posvel,prec_hi,prec_lo);
1808 else if (magic_int==MAGIC_INT_VEL0x56474E54)
1809 return tng_compress_uncompress_vel_int(data,posvel,prec_hi,prec_lo);
1810 else
1811 return 1;
1812}
1813
1814void DECLSPECDLLEXPORT tng_compress_int_to_double(int *posvel_int,unsigned long prec_hi, unsigned long prec_lo,
1815 int natoms,int nframes,
1816 double *posvel_double)
1817{
1818 unquantize(posvel_double,natoms,nframes,PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),posvel_int);
1819}
1820
1821void DECLSPECDLLEXPORT tng_compress_int_to_float(int *posvel_int,unsigned long prec_hi, unsigned long prec_lo,
1822 int natoms,int nframes,
1823 float *posvel_float)
1824{
1825 unquantize_float(posvel_float,natoms,nframes,(float)PRECISION(prec_hi,prec_lo)(Ptngc_i32x2_to_d(prec_hi,prec_lo)),posvel_int);
1826}
1827
1828static char *compress_algo_pos[TNG_COMPRESS_ALGO_MAX11]={
1829 "Positions invalid algorithm",
1830 "Positions stopbits interframe",
1831 "Positions triplet interframe",
1832 "Positions triplet intraframe",
1833 "Positions invalid algorithm",
1834 "Positions XTC2",
1835 "Positions invalid algorithm",
1836 "Positions triplet one to one",
1837 "Positions BWLZH interframe",
1838 "Positions BWLZH intraframe",
1839 "Positions XTC3"
1840};
1841
1842static char *compress_algo_vel[TNG_COMPRESS_ALGO_MAX11]={
1843 "Velocities invalid algorithm",
1844 "Velocities stopbits one to one",
1845 "Velocities triplet interframe",
1846 "Velocities triplet one to one",
1847 "Velocities invalid algorithm",
1848 "Velocities invalid algorithm",
1849 "Velocities stopbits interframe",
1850 "Velocities invalid algorithm",
1851 "Velocities BWLZH interframe",
1852 "Velocities BWLZH one to one",
1853 "Velocities invalid algorithm"
1854};
1855
1856char DECLSPECDLLEXPORT *tng_compress_initial_pos_algo(int *algo)
1857{
1858 int i=algo[0];
1859 if (i<0)
1860 i=0;
1861 if (i>=TNG_COMPRESS_ALGO_MAX11)
1862 i=0;
1863 return compress_algo_pos[i];
1864}
1865
1866char DECLSPECDLLEXPORT *tng_compress_pos_algo(int *algo)
1867{
1868 int i=algo[2];
1869 if (i<0)
1870 i=0;
1871 if (i>=TNG_COMPRESS_ALGO_MAX11)
1872 i=0;
1873 return compress_algo_pos[i];
1874}
1875
1876char DECLSPECDLLEXPORT *tng_compress_initial_vel_algo(int *algo)
1877{
1878 int i=algo[0];
1879 if (i<0)
1880 i=0;
1881 if (i>=TNG_COMPRESS_ALGO_MAX11)
1882 i=0;
1883 return compress_algo_vel[i];
1884}
1885
1886char DECLSPECDLLEXPORT *tng_compress_vel_algo(int *algo)
1887{
1888 int i=algo[2];
1889 if (i<0)
1890 i=0;
1891 if (i>=TNG_COMPRESS_ALGO_MAX11)
1892 i=0;
1893 return compress_algo_vel[i];
1894}