1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: t; c-basic-offset: 4 -*-
2 *
3 * Copyright (c) 2009-2014, Erik Lindahl & David van der Spoel
4 * All rights reserved.
5 *
6 * Redistribution and use in source and binary forms, with or without
7 * modification, are permitted provided that the following conditions are met:
8 *
9 * 1. Redistributions of source code must retain the above copyright notice, this
10 * list of conditions and the following disclaimer.
11 *
12 * 2. Redistributions in binary form must reproduce the above copyright notice,
13 * this list of conditions and the following disclaimer in the documentation
14 * and/or other materials provided with the distribution.
15 *
16 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
19 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
20 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
21 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
22 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
23 * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
24 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 */
27
28 /* Get HAVE_RPC_XDR_H, F77_FUNC from config.h if available */
29 #ifdef HAVE_CONFIG_H
30 #include <config.h>
31 #endif
32
33 #include <stdio.h>
34 #include <stdlib.h>
35 #include <string.h>
36 #include <math.h>
37 #include <limits.h>
38 #include <errno.h>
39
40 /* get fixed-width types if we are using ANSI C99 */
41 #ifdef HAVE_STDINT_H
42 # include <stdint.h>
43 #elif (defined HAVE_INTTYPES_H)
44 # include <inttypes.h>
45 #endif
46
47 #include <sys/types.h>
48 #define _FILE_OFFSET_BITS 64
49
50 #ifdef _WIN32
51 # include <io.h>
52 #else /* __unix__ */
53 # include <unistd.h>
54 #endif
55
56
57 #ifdef HAVE_RPC_XDR_H
58 # include <rpc/rpc.h>
59 # include <rpc/xdr.h>
60 #endif
61
62 #include "xdrfile.h"
63
64 /* Default FORTRAN name mangling is: lower case name, append underscore */
65 #ifndef F77_FUNC
66 #define F77_FUNC(name,NAME) name ## _
67 #endif
68
69 char *exdr_message[exdrNR] = {
70 "OK",
71 "Header",
72 "String",
73 "Double",
74 "Integer",
75 "Float",
76 "Unsigned integer",
77 "Compressed 3D coordinate",
78 "Closing file",
79 "Magic number",
80 "Not enough memory",
81 "End of file",
82 "File not found"
83 };
84
85 /*
86 * Declare our own XDR routines statically if no libraries are present.
87 * Actual implementation is at the end of this file.
88 *
89 * We don't want the low-level XDR implementation as part of the Gromacs
90 * documentation, so skip it for doxygen too...
91 */
92 #if (!defined HAVE_RPC_XDR_H && !defined DOXYGEN)
93
94 enum xdr_op
95 {
96 XDR_ENCODE = 0,
97 XDR_DECODE = 1,
98 XDR_FREE = 2
99 };
100
101
102 /* We need integer types that are guaranteed to be 4 bytes wide.
103 * If ANSI C99 headers were included they are already defined
104 * as int32_t and uint32_t. Check, and if not define them ourselves.
105 * Since it is just our workaround for missing ANSI C99 types, avoid adding
106 * it to the doxygen documentation.
107 */
108 #if !(defined INT32_MAX || defined DOXYGEN)
109 # if (INT_MAX == 2147483647)
110 # define int32_t int
111 # define uint32_t unsigned int
112 # define INT32_MAX 2147483647
113 # elif (LONG_MAX == 2147483647)
114 # define int32_t long
115 # define uint32_t unsigned long
116 # define INT32_MAX 2147483647L
117 # else
118 # error ERROR: No 32 bit wide integer type found!
119 # error Use system XDR libraries instead, or update xdrfile.c
120 # endif
121 #endif
122
123 typedef struct XDR XDR;
124
125 struct XDR
126 {
127 enum xdr_op x_op;
128 struct xdr_ops
129 {
130 int (*x_getlong) (XDR *__xdrs, int32_t *__lp);
131 int (*x_putlong) (XDR *__xdrs, int32_t *__lp);
132 int (*x_getbytes) (XDR *__xdrs, char *__addr, unsigned int __len);
133 int (*x_putbytes) (XDR *__xdrs, char *__addr, unsigned int __len);
134 /* two next routines are not 64-bit IO safe - don't use! */
135 unsigned int (*x_getpostn) (XDR *__xdrs);
136 int (*x_setpostn) (XDR *__xdrs, unsigned int __pos);
137 void (*x_destroy) (XDR *__xdrs);
138 }
139 *x_ops;
140 char *x_private;
141 };
142
143 static int xdr_char (XDR *xdrs, char *ip);
144 static int xdr_u_char (XDR *xdrs, unsigned char *ip);
145 static int xdr_short (XDR *xdrs, short *ip);
146 static int xdr_u_short (XDR *xdrs, unsigned short *ip);
147 static int xdr_int (XDR *xdrs, int *ip);
148 static int xdr_u_int (XDR *xdrs, unsigned int *ip);
149 static int xdr_float (XDR *xdrs, float *ip);
150 static int xdr_double (XDR *xdrs, double *ip);
151 static int xdr_string (XDR *xdrs, char **ip, unsigned int maxsize);
152 static int xdr_opaque (XDR *xdrs, char *cp, unsigned int cnt);
153 static void xdrstdio_create (XDR *xdrs, FILE *fp, enum xdr_op xop);
154
155 #define xdr_getpos(xdrs) \
156 (*(xdrs)->x_ops->x_getpostn)(xdrs)
157 #define xdr_setpos(xdrs, pos) \
158 (*(xdrs)->x_ops->x_setpostn)(xdrs, pos)
159 #define xdr_destroy(xdrs) \
160 do { \
161 if ((xdrs)->x_ops->x_destroy) \
162 (*(xdrs)->x_ops->x_destroy)(xdrs); \
163 } while (0)
164 #endif /* end of our own XDR declarations */
165
166
167
168
169
170 /** Contents of the abstract XDRFILE data structure.
171 *
172 * @internal
173 *
174 * This structure is used to provide an XDR file interface that is
175 * virtual identical to the standard UNIX fopen/fread/fwrite/fclose.
176 */
177 struct XDRFILE
178 {
179 FILE * fp; /**< pointer to standard C library file handle */
180 XDR * xdr; /**< pointer to corresponding XDR handle */
181 char mode; /**< r=read, w=write, a=append */
182 int * buf1; /**< Buffer for internal use */
183 int buf1size; /**< Current allocated length of buf1 */
184 int * buf2; /**< Buffer for internal use */
185 int buf2size; /**< Current allocated length of buf2 */
186 };
187
188
189
190
191 /*************************************************************
192 * Implementation of higher-level routines to read/write *
193 * portable data based on the XDR standard. These should be *
194 * called from C - see further down for Fortran77 wrappers. *
195 *************************************************************/
196
197 XDRFILE *
xdrfile_open(const char * path,const char * mode)198 xdrfile_open(const char *path, const char *mode)
199 {
200 char newmode[5];
201 enum xdr_op xdrmode;
202 XDRFILE *xfp;
203
204 /* make sure XDR files are opened in binary mode... */
205 if(*mode=='w' || *mode=='W')
206 {
207 sprintf(newmode,"wb+");
208 xdrmode=XDR_ENCODE;
209 } else if(*mode == 'a' || *mode == 'A')
210 {
211 sprintf(newmode,"ab+");
212 xdrmode = XDR_ENCODE;
213 } else if(*mode == 'r' || *mode == 'R')
214 {
215 sprintf(newmode,"rb");
216 xdrmode = XDR_DECODE;
217 } else /* cannot determine mode */
218 return NULL;
219
220 if((xfp=(XDRFILE *)malloc(sizeof(XDRFILE)))==NULL)
221 return NULL;
222 if((xfp->fp=fopen(path,newmode))==NULL)
223 {
224 free(xfp);
225 return NULL;
226 }
227 if((xfp->xdr=(XDR *)malloc(sizeof(XDR)))==NULL)
228 {
229 fclose(xfp->fp);
230 free(xfp);
231 return NULL;
232 }
233 xfp->mode=*mode;
234 xdrstdio_create((XDR *)(xfp->xdr),xfp->fp,xdrmode);
235 xfp->buf1 = xfp->buf2 = NULL;
236 xfp->buf1size = xfp->buf2size = 0;
237 return xfp;
238 }
239
240 int
xdrfile_close(XDRFILE * xfp)241 xdrfile_close(XDRFILE *xfp)
242 {
243 int ret=exdrCLOSE;
244 if(xfp)
245 {
246 /* flush and destroy XDR stream */
247 if(xfp->xdr)
248 xdr_destroy((XDR *)(xfp->xdr));
249 free(xfp->xdr);
250 /* close the file */
251 ret=fclose(xfp->fp);
252 if(xfp->buf1size)
253 free(xfp->buf1);
254 if(xfp->buf2size)
255 free(xfp->buf2);
256 free(xfp);
257 }
258 return ret; /* return 0 if ok */
259 }
260
261
262
263 int
xdrfile_read_int(int * ptr,int ndata,XDRFILE * xfp)264 xdrfile_read_int(int *ptr, int ndata, XDRFILE* xfp)
265 {
266 int i=0;
267
268 /* read write is encoded in the XDR struct */
269 while(i<ndata && xdr_int((XDR *)(xfp->xdr),ptr+i))
270 i++;
271
272 return i;
273 }
274
275 int
xdrfile_write_int(int * ptr,int ndata,XDRFILE * xfp)276 xdrfile_write_int(int *ptr, int ndata, XDRFILE* xfp)
277 {
278 int i=0;
279
280 /* read write is encoded in the XDR struct */
281 while(i<ndata && xdr_int((XDR *)(xfp->xdr),ptr+i))
282 i++;
283 return i;
284 }
285
286
287 int
xdrfile_read_uint(unsigned int * ptr,int ndata,XDRFILE * xfp)288 xdrfile_read_uint(unsigned int *ptr, int ndata, XDRFILE* xfp)
289 {
290 int i=0;
291
292 /* read write is encoded in the XDR struct */
293 while(i<ndata && xdr_u_int((XDR *)(xfp->xdr),ptr+i))
294 i++;
295
296 return i;
297 }
298
299 int
xdrfile_write_uint(unsigned int * ptr,int ndata,XDRFILE * xfp)300 xdrfile_write_uint(unsigned int *ptr, int ndata, XDRFILE* xfp)
301 {
302 int i=0;
303
304 /* read write is encoded in the XDR struct */
305 while(i<ndata && xdr_u_int((XDR *)(xfp->xdr),ptr+i))
306 i++;
307 return i;
308 }
309
310 int
xdrfile_read_char(char * ptr,int ndata,XDRFILE * xfp)311 xdrfile_read_char(char *ptr, int ndata, XDRFILE* xfp)
312 {
313 int i=0;
314
315 /* read write is encoded in the XDR struct */
316 while(i<ndata && xdr_char((XDR *)(xfp->xdr),ptr+i))
317 i++;
318
319 return i;
320 }
321
322 int
xdrfile_write_char(char * ptr,int ndata,XDRFILE * xfp)323 xdrfile_write_char(char *ptr, int ndata, XDRFILE* xfp)
324 {
325 int i=0;
326
327 /* read write is encoded in the XDR struct */
328 while(i<ndata && xdr_char((XDR *)(xfp->xdr),ptr+i))
329 i++;
330 return i;
331 }
332
333
334 int
xdrfile_read_uchar(unsigned char * ptr,int ndata,XDRFILE * xfp)335 xdrfile_read_uchar(unsigned char *ptr, int ndata, XDRFILE* xfp)
336 {
337 int i=0;
338
339 /* read write is encoded in the XDR struct */
340 while(i<ndata && xdr_u_char((XDR *)(xfp->xdr),ptr+i))
341 i++;
342
343 return i;
344 }
345
346 int
xdrfile_write_uchar(unsigned char * ptr,int ndata,XDRFILE * xfp)347 xdrfile_write_uchar(unsigned char *ptr, int ndata, XDRFILE* xfp)
348 {
349 int i=0;
350
351 /* read write is encoded in the XDR struct */
352 while(i<ndata && xdr_u_char((XDR *)(xfp->xdr),ptr+i))
353 i++;
354 return i;
355 }
356
357 int
xdrfile_read_short(short * ptr,int ndata,XDRFILE * xfp)358 xdrfile_read_short(short *ptr, int ndata, XDRFILE* xfp)
359 {
360 int i=0;
361
362 /* read write is encoded in the XDR struct */
363 while(i<ndata && xdr_short((XDR *)(xfp->xdr),ptr+i))
364 i++;
365
366 return i;
367 }
368
369 int
xdrfile_write_short(short * ptr,int ndata,XDRFILE * xfp)370 xdrfile_write_short(short *ptr, int ndata, XDRFILE* xfp)
371 {
372 int i=0;
373
374 /* read write is encoded in the XDR struct */
375 while(i<ndata && xdr_short((XDR *)(xfp->xdr),ptr+i))
376 i++;
377 return i;
378 }
379
380
381 int
xdrfile_read_ushort(unsigned short * ptr,int ndata,XDRFILE * xfp)382 xdrfile_read_ushort(unsigned short *ptr, int ndata, XDRFILE* xfp)
383 {
384 int i=0;
385
386 /* read write is encoded in the XDR struct */
387 while(i<ndata && xdr_u_short((XDR *)(xfp->xdr),ptr+i))
388 i++;
389
390 return i;
391 }
392
393 int
xdrfile_write_ushort(unsigned short * ptr,int ndata,XDRFILE * xfp)394 xdrfile_write_ushort(unsigned short *ptr, int ndata, XDRFILE* xfp)
395 {
396 int i=0;
397
398 /* read write is encoded in the XDR struct */
399 while(i<ndata && xdr_u_short((XDR *)(xfp->xdr),ptr+i))
400 i++;
401 return i;
402 }
403
404 int
xdrfile_read_float(float * ptr,int ndata,XDRFILE * xfp)405 xdrfile_read_float(float *ptr, int ndata, XDRFILE* xfp)
406 {
407 int i=0;
408 /* read write is encoded in the XDR struct */
409 while(i<ndata && xdr_float((XDR *)(xfp->xdr),ptr+i))
410 i++;
411 return i;
412 }
413
414 int
xdrfile_write_float(float * ptr,int ndata,XDRFILE * xfp)415 xdrfile_write_float(float *ptr, int ndata, XDRFILE* xfp)
416 {
417 int i=0;
418 /* read write is encoded in the XDR struct */
419 while(i<ndata && xdr_float((XDR *)(xfp->xdr),ptr+i))
420 i++;
421 return i;
422 }
423
424 int
xdrfile_read_double(double * ptr,int ndata,XDRFILE * xfp)425 xdrfile_read_double(double *ptr, int ndata, XDRFILE* xfp)
426 {
427 int i=0;
428 /* read write is encoded in the XDR struct */
429 while(i<ndata && xdr_double((XDR *)(xfp->xdr),ptr+i))
430 i++;
431 return i;
432 }
433
434 int
xdrfile_write_double(double * ptr,int ndata,XDRFILE * xfp)435 xdrfile_write_double(double *ptr, int ndata, XDRFILE* xfp)
436 {
437 int i=0;
438 /* read write is encoded in the XDR struct */
439 while(i<ndata && xdr_double((XDR *)(xfp->xdr),ptr+i))
440 i++;
441 return i;
442 }
443
444 int
xdrfile_read_string(char * ptr,int maxlen,XDRFILE * xfp)445 xdrfile_read_string(char *ptr, int maxlen, XDRFILE* xfp)
446 {
447 int i;
448 if(xdr_string((XDR *)(xfp->xdr),&ptr,maxlen)) {
449 i=0;
450 while(i<maxlen && ptr[i]!=0)
451 i++;
452 if(i==maxlen)
453 return maxlen;
454 else
455 return i+1;
456 } else
457 return 0;
458 }
459
460 int
xdrfile_write_string(char * ptr,XDRFILE * xfp)461 xdrfile_write_string(char *ptr, XDRFILE* xfp)
462 {
463 int len=strlen(ptr)+1;
464
465 if(xdr_string((XDR *)(xfp->xdr),&ptr,len))
466 return len;
467 else
468 return 0;
469 }
470
471
472 int
xdrfile_read_opaque(char * ptr,int cnt,XDRFILE * xfp)473 xdrfile_read_opaque(char *ptr, int cnt, XDRFILE* xfp)
474 {
475 if(xdr_opaque((XDR *)(xfp->xdr),ptr,cnt))
476 return cnt;
477 else
478 return 0;
479 }
480
481
482 int
xdrfile_write_opaque(char * ptr,int cnt,XDRFILE * xfp)483 xdrfile_write_opaque(char *ptr, int cnt, XDRFILE* xfp)
484 {
485 if(xdr_opaque((XDR *)(xfp->xdr),ptr,cnt))
486 return cnt;
487 else
488 return 0;
489 }
490
491
492 /* Internal support routines for reading/writing compressed coordinates
493 * sizeofint - calculate smallest number of bits necessary
494 * to represent a certain integer.
495 */
496 static int
sizeofint(int size)497 sizeofint(int size) {
498 unsigned int num = 1;
499 int num_of_bits = 0;
500
501 while (size >= num && num_of_bits < 32)
502 {
503 num_of_bits++;
504 num <<= 1;
505 }
506 return num_of_bits;
507 }
508
509
510 /*
511 * sizeofints - calculate 'bitsize' of compressed ints
512 *
513 * given a number of small unsigned integers and the maximum value
514 * return the number of bits needed to read or write them with the
515 * routines encodeints/decodeints. You need this parameter when
516 * calling those routines.
517 * (However, in some cases we can just use the variable 'smallidx'
518 * which is the exact number of bits, and them we dont need to call
519 * this routine).
520 */
521 static int
sizeofints(int num_of_ints,unsigned int sizes[])522 sizeofints(int num_of_ints, unsigned int sizes[])
523 {
524 int i, num;
525 unsigned int num_of_bytes, num_of_bits, bytes[32], bytecnt, tmp;
526 num_of_bytes = 1;
527 bytes[0] = 1;
528 num_of_bits = 0;
529 for (i=0; i < num_of_ints; i++)
530 {
531 tmp = 0;
532 for (bytecnt = 0; bytecnt < num_of_bytes; bytecnt++)
533 {
534 tmp = bytes[bytecnt] * sizes[i] + tmp;
535 bytes[bytecnt] = tmp & 0xff;
536 tmp >>= 8;
537 }
538 while (tmp != 0)
539 {
540 bytes[bytecnt++] = tmp & 0xff;
541 tmp >>= 8;
542 }
543 num_of_bytes = bytecnt;
544 }
545 num = 1;
546 num_of_bytes--;
547 while (bytes[num_of_bytes] >= num)
548 {
549 num_of_bits++;
550 num *= 2;
551 }
552 return num_of_bits + num_of_bytes * 8;
553
554 }
555
556
557 /*
558 * encodebits - encode num into buf using the specified number of bits
559 *
560 * This routines appends the value of num to the bits already present in
561 * the array buf. You need to give it the number of bits to use and you had
562 * better make sure that this number of bits is enough to hold the value.
563 * Num must also be positive.
564 */
565 static void
encodebits(int buf[],int num_of_bits,int num)566 encodebits(int buf[], int num_of_bits, int num)
567 {
568
569 unsigned int cnt, lastbyte;
570 int lastbits;
571 unsigned char * cbuf;
572
573 cbuf = ((unsigned char *)buf) + 3 * sizeof(*buf);
574 cnt = (unsigned int) buf[0];
575 lastbits = buf[1];
576 lastbyte =(unsigned int) buf[2];
577 while (num_of_bits >= 8)
578 {
579 lastbyte = (lastbyte << 8) | ((num >> (num_of_bits -8)) /* & 0xff*/);
580 cbuf[cnt++] = lastbyte >> lastbits;
581 num_of_bits -= 8;
582 }
583 if (num_of_bits > 0)
584 {
585 lastbyte = (lastbyte << num_of_bits) | num;
586 lastbits += num_of_bits;
587 if (lastbits >= 8)
588 {
589 lastbits -= 8;
590 cbuf[cnt++] = lastbyte >> lastbits;
591 }
592 }
593 buf[0] = cnt;
594 buf[1] = lastbits;
595 buf[2] = lastbyte;
596 if (lastbits>0)
597 {
598 cbuf[cnt] = lastbyte << (8 - lastbits);
599 }
600 }
601
602 /*
603 * encodeints - encode a small set of small integers in compressed format
604 *
605 * this routine is used internally by xdr3dfcoord, to encode a set of
606 * small integers to the buffer for writing to a file.
607 * Multiplication with fixed (specified maximum) sizes is used to get
608 * to one big, multibyte integer. Allthough the routine could be
609 * modified to handle sizes bigger than 16777216, or more than just
610 * a few integers, this is not done because the gain in compression
611 * isn't worth the effort. Note that overflowing the multiplication
612 * or the byte buffer (32 bytes) is unchecked and whould cause bad results.
613 * THese things are checked in the calling routines, so make sure not
614 * to remove those checks...
615 */
616
617 static void
encodeints(int buf[],int num_of_ints,int num_of_bits,unsigned int sizes[],unsigned int nums[])618 encodeints(int buf[], int num_of_ints, int num_of_bits,
619 unsigned int sizes[], unsigned int nums[])
620 {
621
622 int i;
623 unsigned int bytes[32], num_of_bytes, bytecnt, tmp;
624
625 tmp = nums[0];
626 num_of_bytes = 0;
627 do
628 {
629 bytes[num_of_bytes++] = tmp & 0xff;
630 tmp >>= 8;
631 } while (tmp != 0);
632
633 for (i = 1; i < num_of_ints; i++)
634 {
635 if (nums[i] >= sizes[i])
636 {
637 fprintf(stderr,"major breakdown in encodeints - num %u doesn't "
638 "match size %u\n", nums[i], sizes[i]);
639 abort();
640 }
641 /* use one step multiply */
642 tmp = nums[i];
643 for (bytecnt = 0; bytecnt < num_of_bytes; bytecnt++)
644 {
645 tmp = bytes[bytecnt] * sizes[i] + tmp;
646 bytes[bytecnt] = tmp & 0xff;
647 tmp >>= 8;
648 }
649 while (tmp != 0)
650 {
651 bytes[bytecnt++] = tmp & 0xff;
652 tmp >>= 8;
653 }
654 num_of_bytes = bytecnt;
655 }
656 if (num_of_bits >= num_of_bytes * 8)
657 {
658 for (i = 0; i < num_of_bytes; i++)
659 {
660 encodebits(buf, 8, bytes[i]);
661 }
662 encodebits(buf, num_of_bits - num_of_bytes * 8, 0);
663 }
664 else
665 {
666 for (i = 0; i < num_of_bytes-1; i++)
667 {
668 encodebits(buf, 8, bytes[i]);
669 }
670 encodebits(buf, num_of_bits- (num_of_bytes -1) * 8, bytes[i]);
671 }
672 }
673
674
675 /*
676 * decodebits - decode number from buf using specified number of bits
677 *
678 * extract the number of bits from the array buf and construct an integer
679 * from it. Return that value.
680 *
681 */
682
683 static int
decodebits(int buf[],int num_of_bits)684 decodebits(int buf[], int num_of_bits)
685 {
686
687 int cnt, num;
688 unsigned int lastbits, lastbyte;
689 unsigned char * cbuf;
690 int mask = (1 << num_of_bits) -1;
691
692 cbuf = ((unsigned char *)buf) + 3 * sizeof(*buf);
693 cnt = buf[0];
694 lastbits = (unsigned int) buf[1];
695 lastbyte = (unsigned int) buf[2];
696
697 num = 0;
698 while (num_of_bits >= 8)
699 {
700 lastbyte = ( lastbyte << 8 ) | cbuf[cnt++];
701 num |= (lastbyte >> lastbits) << (num_of_bits - 8);
702 num_of_bits -=8;
703 }
704 if (num_of_bits > 0)
705 {
706 if (lastbits < num_of_bits)
707 {
708 lastbits += 8;
709 lastbyte = (lastbyte << 8) | cbuf[cnt++];
710 }
711 lastbits -= num_of_bits;
712 num |= (lastbyte >> lastbits) & ((1 << num_of_bits) -1);
713 }
714 num &= mask;
715 buf[0] = cnt;
716 buf[1] = lastbits;
717 buf[2] = lastbyte;
718 return num;
719 }
720
721 /*
722 * decodeints - decode 'small' integers from the buf array
723 *
724 * this routine is the inverse from encodeints() and decodes the small integers
725 * written to buf by calculating the remainder and doing divisions with
726 * the given sizes[]. You need to specify the total number of bits to be
727 * used from buf in num_of_bits.
728 *
729 */
730
731 static void
decodeints(int buf[],int num_of_ints,int num_of_bits,unsigned int sizes[],int nums[])732 decodeints(int buf[], int num_of_ints, int num_of_bits,
733 unsigned int sizes[], int nums[])
734 {
735
736 int bytes[32];
737 int i, j, num_of_bytes, p, num;
738
739 bytes[1] = bytes[2] = bytes[3] = 0;
740 num_of_bytes = 0;
741 while (num_of_bits > 8)
742 {
743 bytes[num_of_bytes++] = decodebits(buf, 8);
744 num_of_bits -= 8;
745 }
746 if (num_of_bits > 0)
747 {
748 bytes[num_of_bytes++] = decodebits(buf, num_of_bits);
749 }
750 for (i = num_of_ints-1; i > 0; i--)
751 {
752 num = 0;
753 for (j = num_of_bytes-1; j >=0; j--)
754 {
755 num = (num << 8) | bytes[j];
756 p = num / sizes[i];
757 bytes[j] = p;
758 num = num - p * sizes[i];
759 }
760 nums[i] = num;
761 }
762 nums[0] = bytes[0] | (bytes[1] << 8) | (bytes[2] << 16) | (bytes[3] << 24);
763 }
764
765
766 static const int magicints[] =
767 {
768 0, 0, 0, 0, 0, 0, 0, 0, 0, 8, 10, 12, 16, 20, 25, 32, 40, 50, 64,
769 80, 101, 128, 161, 203, 256, 322, 406, 512, 645, 812, 1024, 1290,
770 1625, 2048, 2580, 3250, 4096, 5060, 6501, 8192, 10321, 13003,
771 16384, 20642, 26007, 32768, 41285, 52015, 65536,82570, 104031,
772 131072, 165140, 208063, 262144, 330280, 416127, 524287, 660561,
773 832255, 1048576, 1321122, 1664510, 2097152, 2642245, 3329021,
774 4194304, 5284491, 6658042, 8388607, 10568983, 13316085, 16777216
775 };
776
777 #define FIRSTIDX 9
778 /* note that magicints[FIRSTIDX-1] == 0 */
779 #define LASTIDX (sizeof(magicints) / sizeof(*magicints))
780
781 /* Compressed coordinate routines - modified from the original
782 * implementation by Frans v. Hoesel to make them threadsafe.
783 */
784 int
xdrfile_decompress_coord_float(float * ptr,int * size,float * precision,XDRFILE * xfp)785 xdrfile_decompress_coord_float(float *ptr,
786 int *size,
787 float *precision,
788 XDRFILE* xfp)
789 {
790 int minint[3], maxint[3], *lip;
791 int smallidx, minidx, maxidx;
792 unsigned sizeint[3], sizesmall[3], bitsizeint[3], size3;
793 int k, *buf1, *buf2, lsize, flag;
794 int smallnum, smaller, larger, i, is_smaller, run;
795 float *lfp, inv_precision;
796 int tmp, *thiscoord, prevcoord[3];
797 unsigned int bitsize;
798
799 bitsizeint[0] = 0;
800 bitsizeint[1] = 0;
801 bitsizeint[2] = 0;
802
803 if(xfp==NULL || ptr==NULL)
804 return -1;
805 tmp=xdrfile_read_int(&lsize,1,xfp);
806 if(tmp==0)
807 return -1; /* return if we could not read size */
808 if (*size < lsize)
809 {
810 fprintf(stderr, "Requested to decompress %d coords, file contains %d\n",
811 *size, lsize);
812 return -1;
813 }
814 *size = lsize;
815 size3 = *size * 3;
816 if(size3>xfp->buf1size)
817 {
818 if((xfp->buf1=(int *)malloc(sizeof(int)*size3))==NULL)
819 {
820 fprintf(stderr,"Cannot allocate memory for decompressing coordinates.\n");
821 return -1;
822 }
823 xfp->buf1size=size3;
824 xfp->buf2size=size3*1.2;
825 if((xfp->buf2=(int *)malloc(sizeof(int)*xfp->buf2size))==NULL)
826 {
827 fprintf(stderr,"Cannot allocate memory for decompressing coordinates.\n");
828 return -1;
829 }
830 }
831 /* Dont bother with compression for three atoms or less */
832 if(*size<=9)
833 {
834 return xdrfile_read_float(ptr,size3,xfp)/3;
835 /* return number of coords, not floats */
836 }
837 /* Compression-time if we got here. Read precision first */
838 xdrfile_read_float(precision,1,xfp);
839
840 /* avoid repeated pointer dereferencing. */
841 buf1=xfp->buf1;
842 buf2=xfp->buf2;
843 /* buf2[0-2] are special and do not contain actual data */
844 buf2[0] = buf2[1] = buf2[2] = 0;
845 xdrfile_read_int(minint,3,xfp);
846 xdrfile_read_int(maxint,3,xfp);
847
848 sizeint[0] = maxint[0] - minint[0]+1;
849 sizeint[1] = maxint[1] - minint[1]+1;
850 sizeint[2] = maxint[2] - minint[2]+1;
851
852 /* check if one of the sizes is to big to be multiplied */
853 if ((sizeint[0] | sizeint[1] | sizeint[2] ) > 0xffffff)
854 {
855 bitsizeint[0] = sizeofint(sizeint[0]);
856 bitsizeint[1] = sizeofint(sizeint[1]);
857 bitsizeint[2] = sizeofint(sizeint[2]);
858 bitsize = 0; /* flag the use of large sizes */
859 }
860 else
861 {
862 bitsize = sizeofints(3, sizeint);
863 }
864
865 if (xdrfile_read_int(&smallidx,1,xfp) == 0)
866 return 0; /* not sure what has happened here or why we return... */
867 tmp=smallidx+8;
868 maxidx = (LASTIDX<tmp) ? LASTIDX : tmp;
869 minidx = maxidx - 8; /* often this equal smallidx */
870 tmp = smallidx-1;
871 tmp = (FIRSTIDX>tmp) ? FIRSTIDX : tmp;
872 smaller = magicints[tmp] / 2;
873 smallnum = magicints[smallidx] / 2;
874 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx] ;
875 larger = magicints[maxidx];
876
877 /* buf2[0] holds the length in bytes */
878
879 if (xdrfile_read_int(buf2,1,xfp) == 0)
880 return 0;
881 if (xdrfile_read_opaque((char *)&(buf2[3]),(unsigned int)buf2[0],xfp) == 0)
882 return 0;
883 buf2[0] = buf2[1] = buf2[2] = 0;
884
885 lfp = ptr;
886 inv_precision = 1.0 / * precision;
887 run = 0;
888 i = 0;
889 lip = buf1;
890 while ( i < lsize )
891 {
892 thiscoord = (int *)(lip) + i * 3;
893
894 if (bitsize == 0)
895 {
896 thiscoord[0] = decodebits(buf2, bitsizeint[0]);
897 thiscoord[1] = decodebits(buf2, bitsizeint[1]);
898 thiscoord[2] = decodebits(buf2, bitsizeint[2]);
899 }
900 else
901 {
902 decodeints(buf2, 3, bitsize, sizeint, thiscoord);
903 }
904
905 i++;
906 thiscoord[0] += minint[0];
907 thiscoord[1] += minint[1];
908 thiscoord[2] += minint[2];
909
910 prevcoord[0] = thiscoord[0];
911 prevcoord[1] = thiscoord[1];
912 prevcoord[2] = thiscoord[2];
913
914 flag = decodebits(buf2, 1);
915 is_smaller = 0;
916 if (flag == 1)
917 {
918 run = decodebits(buf2, 5);
919 is_smaller = run % 3;
920 run -= is_smaller;
921 is_smaller--;
922 }
923 if (run > 0)
924 {
925 thiscoord += 3;
926 for (k = 0; k < run; k+=3)
927 {
928 decodeints(buf2, 3, smallidx, sizesmall, thiscoord);
929 i++;
930 thiscoord[0] += prevcoord[0] - smallnum;
931 thiscoord[1] += prevcoord[1] - smallnum;
932 thiscoord[2] += prevcoord[2] - smallnum;
933 if (k == 0) {
934 /* interchange first with second atom for better
935 * compression of water molecules
936 */
937 tmp = thiscoord[0]; thiscoord[0] = prevcoord[0];
938 prevcoord[0] = tmp;
939 tmp = thiscoord[1]; thiscoord[1] = prevcoord[1];
940 prevcoord[1] = tmp;
941 tmp = thiscoord[2]; thiscoord[2] = prevcoord[2];
942 prevcoord[2] = tmp;
943 *lfp++ = prevcoord[0] * inv_precision;
944 *lfp++ = prevcoord[1] * inv_precision;
945 *lfp++ = prevcoord[2] * inv_precision;
946 } else {
947 prevcoord[0] = thiscoord[0];
948 prevcoord[1] = thiscoord[1];
949 prevcoord[2] = thiscoord[2];
950 }
951 *lfp++ = thiscoord[0] * inv_precision;
952 *lfp++ = thiscoord[1] * inv_precision;
953 *lfp++ = thiscoord[2] * inv_precision;
954 }
955 }
956 else
957 {
958 *lfp++ = thiscoord[0] * inv_precision;
959 *lfp++ = thiscoord[1] * inv_precision;
960 *lfp++ = thiscoord[2] * inv_precision;
961 }
962 smallidx += is_smaller;
963 if (is_smaller < 0)
964 {
965 smallnum = smaller;
966
967 if (smallidx > FIRSTIDX)
968 {
969 smaller = magicints[smallidx - 1] /2;
970 }
971 else
972 {
973 smaller = 0;
974 }
975 }
976 else if (is_smaller > 0)
977 {
978 smaller = smallnum;
979 smallnum = magicints[smallidx] / 2;
980 }
981 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx] ;
982 }
983 return *size;
984 }
985
986 int
xdrfile_compress_coord_float(float * ptr,int size,float precision,XDRFILE * xfp)987 xdrfile_compress_coord_float(float *ptr,
988 int size,
989 float precision,
990 XDRFILE* xfp)
991 {
992 int minint[3], maxint[3], mindiff, *lip, diff;
993 int lint1, lint2, lint3, oldlint1, oldlint2, oldlint3, smallidx;
994 int minidx, maxidx;
995 unsigned sizeint[3], sizesmall[3], bitsizeint[3], size3, *luip;
996 int k, *buf1, *buf2;
997 int smallnum, smaller, larger, i, j, is_small, is_smaller, run, prevrun;
998 float *lfp, lf;
999 int tmp, tmpsum, *thiscoord, prevcoord[3];
1000 unsigned int tmpcoord[30];
1001 int errval=1;
1002 unsigned int bitsize;
1003
1004 if(xfp==NULL)
1005 return -1;
1006 size3=3*size;
1007
1008 bitsizeint[0] = 0;
1009 bitsizeint[1] = 0;
1010 bitsizeint[2] = 0;
1011
1012 if(size3>xfp->buf1size)
1013 {
1014 if((xfp->buf1=(int *)malloc(sizeof(int)*size3))==NULL)
1015 {
1016 fprintf(stderr,"Cannot allocate memory for compressing coordinates.\n");
1017 return -1;
1018 }
1019 xfp->buf1size=size3;
1020 xfp->buf2size=size3*1.2;
1021 if((xfp->buf2=(int *)malloc(sizeof(int)*xfp->buf2size))==NULL)
1022 {
1023 fprintf(stderr,"Cannot allocate memory for compressing coordinates.\n");
1024 return -1;
1025 }
1026 }
1027 if(xdrfile_write_int(&size,1,xfp)==0)
1028 return -1; /* return if we could not write size */
1029 /* Dont bother with compression for three atoms or less */
1030 if(size<=9)
1031 {
1032 return xdrfile_write_float(ptr,size3,xfp)/3;
1033 /* return number of coords, not floats */
1034 }
1035 /* Compression-time if we got here. Write precision first */
1036 if (precision <= 0)
1037 precision = 1000;
1038 xdrfile_write_float(&precision,1,xfp);
1039 /* avoid repeated pointer dereferencing. */
1040 buf1=xfp->buf1;
1041 buf2=xfp->buf2;
1042 /* buf2[0-2] are special and do not contain actual data */
1043 buf2[0] = buf2[1] = buf2[2] = 0;
1044 minint[0] = minint[1] = minint[2] = INT_MAX;
1045 maxint[0] = maxint[1] = maxint[2] = INT_MIN;
1046 prevrun = -1;
1047 lfp = ptr;
1048 lip = buf1;
1049 mindiff = INT_MAX;
1050 oldlint1 = oldlint2 = oldlint3 = 0;
1051 while(lfp < ptr + size3 )
1052 {
1053 /* find nearest integer */
1054 if (*lfp >= 0.0)
1055 lf = *lfp * precision + 0.5;
1056 else
1057 lf = *lfp * precision - 0.5;
1058 if (fabs(lf) > INT_MAX-2)
1059 {
1060 /* scaling would cause overflow */
1061 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1062 errval=0;
1063 }
1064 lint1 = lf;
1065 if (lint1 < minint[0]) minint[0] = lint1;
1066 if (lint1 > maxint[0]) maxint[0] = lint1;
1067 *lip++ = lint1;
1068 lfp++;
1069 if (*lfp >= 0.0)
1070 lf = *lfp * precision + 0.5;
1071 else
1072 lf = *lfp * precision - 0.5;
1073 if (fabs(lf) > INT_MAX-2)
1074 {
1075 /* scaling would cause overflow */
1076 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1077 errval=0;
1078 }
1079 lint2 = lf;
1080 if (lint2 < minint[1]) minint[1] = lint2;
1081 if (lint2 > maxint[1]) maxint[1] = lint2;
1082 *lip++ = lint2;
1083 lfp++;
1084 if (*lfp >= 0.0)
1085 lf = *lfp * precision + 0.5;
1086 else
1087 lf = *lfp * precision - 0.5;
1088 if (fabs(lf) > INT_MAX-2)
1089 {
1090 errval=0;
1091 }
1092 lint3 = lf;
1093 if (lint3 < minint[2]) minint[2] = lint3;
1094 if (lint3 > maxint[2]) maxint[2] = lint3;
1095 *lip++ = lint3;
1096 lfp++;
1097 diff = abs(oldlint1-lint1)+abs(oldlint2-lint2)+abs(oldlint3-lint3);
1098 if (diff < mindiff && lfp > ptr + 3)
1099 mindiff = diff;
1100 oldlint1 = lint1;
1101 oldlint2 = lint2;
1102 oldlint3 = lint3;
1103 }
1104 xdrfile_write_int(minint,3,xfp);
1105 xdrfile_write_int(maxint,3,xfp);
1106
1107 if ((float)maxint[0] - (float)minint[0] >= INT_MAX-2 ||
1108 (float)maxint[1] - (float)minint[1] >= INT_MAX-2 ||
1109 (float)maxint[2] - (float)minint[2] >= INT_MAX-2) {
1110 /* turning value in unsigned by subtracting minint
1111 * would cause overflow
1112 */
1113 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1114 errval=0;
1115 }
1116 sizeint[0] = maxint[0] - minint[0]+1;
1117 sizeint[1] = maxint[1] - minint[1]+1;
1118 sizeint[2] = maxint[2] - minint[2]+1;
1119
1120 /* check if one of the sizes is to big to be multiplied */
1121 if ((sizeint[0] | sizeint[1] | sizeint[2] ) > 0xffffff)
1122 {
1123 bitsizeint[0] = sizeofint(sizeint[0]);
1124 bitsizeint[1] = sizeofint(sizeint[1]);
1125 bitsizeint[2] = sizeofint(sizeint[2]);
1126 bitsize = 0; /* flag the use of large sizes */
1127 }
1128 else
1129 {
1130 bitsize = sizeofints(3, sizeint);
1131 }
1132 lip = buf1;
1133 luip = (unsigned int *) buf1;
1134 smallidx = FIRSTIDX;
1135 while (smallidx < LASTIDX && magicints[smallidx] < mindiff)
1136 {
1137 smallidx++;
1138 }
1139 xdrfile_write_int(&smallidx,1,xfp);
1140 tmp=smallidx+8;
1141 maxidx = (LASTIDX<tmp) ? LASTIDX : tmp;
1142 minidx = maxidx - 8; /* often this equal smallidx */
1143 tmp=smallidx-1;
1144 tmp= (FIRSTIDX>tmp) ? FIRSTIDX : tmp;
1145 smaller = magicints[tmp] / 2;
1146 smallnum = magicints[smallidx] / 2;
1147 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx];
1148 larger = magicints[maxidx] / 2;
1149 i = 0;
1150 while (i < size)
1151 {
1152 is_small = 0;
1153 thiscoord = (int *)(luip) + i * 3;
1154 if (smallidx < maxidx && i >= 1 &&
1155 abs(thiscoord[0] - prevcoord[0]) < larger &&
1156 abs(thiscoord[1] - prevcoord[1]) < larger &&
1157 abs(thiscoord[2] - prevcoord[2]) < larger) {
1158 is_smaller = 1;
1159 }
1160 else if (smallidx > minidx)
1161 {
1162 is_smaller = -1;
1163 }
1164 else
1165 {
1166 is_smaller = 0;
1167 }
1168 if (i + 1 < size)
1169 {
1170 if (abs(thiscoord[0] - thiscoord[3]) < smallnum &&
1171 abs(thiscoord[1] - thiscoord[4]) < smallnum &&
1172 abs(thiscoord[2] - thiscoord[5]) < smallnum)
1173 {
1174 /* interchange first with second atom for better
1175 * compression of water molecules
1176 */
1177 tmp = thiscoord[0]; thiscoord[0] = thiscoord[3];
1178 thiscoord[3] = tmp;
1179 tmp = thiscoord[1]; thiscoord[1] = thiscoord[4];
1180 thiscoord[4] = tmp;
1181 tmp = thiscoord[2]; thiscoord[2] = thiscoord[5];
1182 thiscoord[5] = tmp;
1183 is_small = 1;
1184 }
1185 }
1186 tmpcoord[0] = thiscoord[0] - minint[0];
1187 tmpcoord[1] = thiscoord[1] - minint[1];
1188 tmpcoord[2] = thiscoord[2] - minint[2];
1189 if (bitsize == 0)
1190 {
1191 encodebits(buf2, bitsizeint[0], tmpcoord[0]);
1192 encodebits(buf2, bitsizeint[1], tmpcoord[1]);
1193 encodebits(buf2, bitsizeint[2], tmpcoord[2]);
1194 }
1195 else
1196 {
1197 encodeints(buf2, 3, bitsize, sizeint, tmpcoord);
1198 }
1199 prevcoord[0] = thiscoord[0];
1200 prevcoord[1] = thiscoord[1];
1201 prevcoord[2] = thiscoord[2];
1202 thiscoord = thiscoord + 3;
1203 i++;
1204
1205 run = 0;
1206 if (is_small == 0 && is_smaller == -1)
1207 is_smaller = 0;
1208 while (is_small && run < 8*3)
1209 {
1210 tmpsum=0;
1211 for(j=0;j<3;j++)
1212 {
1213 tmp=thiscoord[j] - prevcoord[j];
1214 tmpsum+=tmp*tmp;
1215 }
1216 if (is_smaller == -1 && tmpsum >= smaller * smaller)
1217 {
1218 is_smaller = 0;
1219 }
1220
1221 tmpcoord[run++] = thiscoord[0] - prevcoord[0] + smallnum;
1222 tmpcoord[run++] = thiscoord[1] - prevcoord[1] + smallnum;
1223 tmpcoord[run++] = thiscoord[2] - prevcoord[2] + smallnum;
1224
1225 prevcoord[0] = thiscoord[0];
1226 prevcoord[1] = thiscoord[1];
1227 prevcoord[2] = thiscoord[2];
1228
1229 i++;
1230 thiscoord = thiscoord + 3;
1231 is_small = 0;
1232 if (i < size &&
1233 abs(thiscoord[0] - prevcoord[0]) < smallnum &&
1234 abs(thiscoord[1] - prevcoord[1]) < smallnum &&
1235 abs(thiscoord[2] - prevcoord[2]) < smallnum)
1236 {
1237 is_small = 1;
1238 }
1239 }
1240 if (run != prevrun || is_smaller != 0)
1241 {
1242 prevrun = run;
1243 encodebits(buf2, 1, 1); /* flag the change in run-length */
1244 encodebits(buf2, 5, run+is_smaller+1);
1245 }
1246 else
1247 {
1248 encodebits(buf2, 1, 0); /* flag the fact that runlength did not change */
1249 }
1250 for (k=0; k < run; k+=3)
1251 {
1252 encodeints(buf2, 3, smallidx, sizesmall, &tmpcoord[k]);
1253 }
1254 if (is_smaller != 0)
1255 {
1256 smallidx += is_smaller;
1257 if (is_smaller < 0)
1258 {
1259 smallnum = smaller;
1260 smaller = magicints[smallidx-1] / 2;
1261 }
1262 else
1263 {
1264 smaller = smallnum;
1265 smallnum = magicints[smallidx] / 2;
1266 }
1267 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx];
1268 }
1269 }
1270 if (buf2[1] != 0) buf2[0]++;
1271 xdrfile_write_int(buf2,1,xfp); /* buf2[0] holds the length in bytes */
1272 tmp=xdrfile_write_opaque((char *)&(buf2[3]),(unsigned int)buf2[0],xfp);
1273 if(tmp==(unsigned int)buf2[0])
1274 return size;
1275 else
1276 return -1;
1277 }
1278
1279
1280 int
xdrfile_decompress_coord_double(double * ptr,int * size,double * precision,XDRFILE * xfp)1281 xdrfile_decompress_coord_double(double *ptr,
1282 int *size,
1283 double *precision,
1284 XDRFILE* xfp)
1285 {
1286 int minint[3], maxint[3], *lip;
1287 int smallidx, minidx, maxidx;
1288 unsigned sizeint[3], sizesmall[3], bitsizeint[3], size3;
1289 int k, *buf1, *buf2, lsize, flag;
1290 int smallnum, smaller, larger, i, is_smaller, run;
1291 double *lfp, inv_precision;
1292 float float_prec, tmpdata[30];
1293 int tmp, *thiscoord, prevcoord[3];
1294 unsigned int bitsize;
1295
1296 bitsizeint[0] = 0;
1297 bitsizeint[1] = 0;
1298 bitsizeint[2] = 0;
1299
1300 if(xfp==NULL || ptr==NULL)
1301 return -1;
1302 tmp=xdrfile_read_int(&lsize,1,xfp);
1303 if(tmp==0)
1304 return -1; /* return if we could not read size */
1305 if (*size < lsize)
1306 {
1307 fprintf(stderr, "Requested to decompress %d coords, file contains %d\n",
1308 *size, lsize);
1309 return -1;
1310 }
1311 *size = lsize;
1312 size3 = *size * 3;
1313 if(size3>xfp->buf1size)
1314 {
1315 if((xfp->buf1=(int *)malloc(sizeof(int)*size3))==NULL)
1316 {
1317 fprintf(stderr,"Cannot allocate memory for decompression coordinates.\n");
1318 return -1;
1319 }
1320 xfp->buf1size=size3;
1321 xfp->buf2size=size3*1.2;
1322 if((xfp->buf2=(int *)malloc(sizeof(int)*xfp->buf2size))==NULL)
1323 {
1324 fprintf(stderr,"Cannot allocate memory for decompressing coordinates.\n");
1325 return -1;
1326 }
1327 }
1328 /* Dont bother with compression for three atoms or less */
1329 if(*size<=9)
1330 {
1331 tmp=xdrfile_read_float(tmpdata,size3,xfp);
1332 for(i=0;i<9*3;i++)
1333 ptr[i]=tmpdata[i];
1334 return tmp/3;
1335 /* return number of coords, not floats */
1336 }
1337 /* Compression-time if we got here. Read precision first */
1338 xdrfile_read_float(&float_prec,1,xfp);
1339 *precision=float_prec;
1340 /* avoid repeated pointer dereferencing. */
1341 buf1=xfp->buf1;
1342 buf2=xfp->buf2;
1343 /* buf2[0-2] are special and do not contain actual data */
1344 buf2[0] = buf2[1] = buf2[2] = 0;
1345 xdrfile_read_int(minint,3,xfp);
1346 xdrfile_read_int(maxint,3,xfp);
1347
1348 sizeint[0] = maxint[0] - minint[0]+1;
1349 sizeint[1] = maxint[1] - minint[1]+1;
1350 sizeint[2] = maxint[2] - minint[2]+1;
1351
1352 /* check if one of the sizes is to big to be multiplied */
1353 if ((sizeint[0] | sizeint[1] | sizeint[2] ) > 0xffffff)
1354 {
1355 bitsizeint[0] = sizeofint(sizeint[0]);
1356 bitsizeint[1] = sizeofint(sizeint[1]);
1357 bitsizeint[2] = sizeofint(sizeint[2]);
1358 bitsize = 0; /* flag the use of large sizes */
1359 }
1360 else
1361 {
1362 bitsize = sizeofints(3, sizeint);
1363 }
1364
1365 if (xdrfile_read_int(&smallidx,1,xfp) == 0)
1366 return 0;
1367 tmp=smallidx+8;
1368 maxidx = (LASTIDX<tmp) ? LASTIDX : tmp;
1369 minidx = maxidx - 8; /* often this equal smallidx */
1370 tmp = smallidx-1;
1371 tmp = (FIRSTIDX>tmp) ? FIRSTIDX : tmp;
1372 smaller = magicints[tmp] / 2;
1373 smallnum = magicints[smallidx] / 2;
1374 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx] ;
1375 larger = magicints[maxidx];
1376
1377 /* buf2[0] holds the length in bytes */
1378
1379 if (xdrfile_read_int(buf2,1,xfp) == 0)
1380 return 0;
1381 if (xdrfile_read_opaque((char *)&(buf2[3]),(unsigned int)buf2[0],xfp) == 0)
1382 return 0;
1383 buf2[0] = buf2[1] = buf2[2] = 0;
1384
1385 lfp = ptr;
1386 inv_precision = 1.0 / * precision;
1387 run = 0;
1388 i = 0;
1389 lip = buf1;
1390 while ( i < lsize )
1391 {
1392 thiscoord = (int *)(lip) + i * 3;
1393
1394 if (bitsize == 0)
1395 {
1396 thiscoord[0] = decodebits(buf2, bitsizeint[0]);
1397 thiscoord[1] = decodebits(buf2, bitsizeint[1]);
1398 thiscoord[2] = decodebits(buf2, bitsizeint[2]);
1399 } else {
1400 decodeints(buf2, 3, bitsize, sizeint, thiscoord);
1401 }
1402
1403 i++;
1404 thiscoord[0] += minint[0];
1405 thiscoord[1] += minint[1];
1406 thiscoord[2] += minint[2];
1407
1408 prevcoord[0] = thiscoord[0];
1409 prevcoord[1] = thiscoord[1];
1410 prevcoord[2] = thiscoord[2];
1411
1412 flag = decodebits(buf2, 1);
1413 is_smaller = 0;
1414 if (flag == 1)
1415 {
1416 run = decodebits(buf2, 5);
1417 is_smaller = run % 3;
1418 run -= is_smaller;
1419 is_smaller--;
1420 }
1421 if (run > 0)
1422 {
1423 thiscoord += 3;
1424 for (k = 0; k < run; k+=3)
1425 {
1426 decodeints(buf2, 3, smallidx, sizesmall, thiscoord);
1427 i++;
1428 thiscoord[0] += prevcoord[0] - smallnum;
1429 thiscoord[1] += prevcoord[1] - smallnum;
1430 thiscoord[2] += prevcoord[2] - smallnum;
1431 if (k == 0)
1432 {
1433 /* interchange first with second atom for better
1434 * compression of water molecules
1435 */
1436 tmp = thiscoord[0]; thiscoord[0] = prevcoord[0];
1437 prevcoord[0] = tmp;
1438 tmp = thiscoord[1]; thiscoord[1] = prevcoord[1];
1439 prevcoord[1] = tmp;
1440 tmp = thiscoord[2]; thiscoord[2] = prevcoord[2];
1441 prevcoord[2] = tmp;
1442 *lfp++ = prevcoord[0] * inv_precision;
1443 *lfp++ = prevcoord[1] * inv_precision;
1444 *lfp++ = prevcoord[2] * inv_precision;
1445 }
1446 else
1447 {
1448 prevcoord[0] = thiscoord[0];
1449 prevcoord[1] = thiscoord[1];
1450 prevcoord[2] = thiscoord[2];
1451 }
1452 *lfp++ = thiscoord[0] * inv_precision;
1453 *lfp++ = thiscoord[1] * inv_precision;
1454 *lfp++ = thiscoord[2] * inv_precision;
1455 }
1456 } else {
1457 *lfp++ = thiscoord[0] * inv_precision;
1458 *lfp++ = thiscoord[1] * inv_precision;
1459 *lfp++ = thiscoord[2] * inv_precision;
1460 }
1461 smallidx += is_smaller;
1462 if (is_smaller < 0) {
1463 smallnum = smaller;
1464 if (smallidx > FIRSTIDX) {
1465 smaller = magicints[smallidx - 1] /2;
1466 } else {
1467 smaller = 0;
1468 }
1469 } else if (is_smaller > 0) {
1470 smaller = smallnum;
1471 smallnum = magicints[smallidx] / 2;
1472 }
1473 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx] ;
1474 }
1475 return *size;
1476 }
1477
1478 int
xdrfile_compress_coord_double(double * ptr,int size,double precision,XDRFILE * xfp)1479 xdrfile_compress_coord_double(double *ptr,
1480 int size,
1481 double precision,
1482 XDRFILE* xfp)
1483 {
1484 int minint[3], maxint[3], mindiff, *lip, diff;
1485 int lint1, lint2, lint3, oldlint1, oldlint2, oldlint3, smallidx;
1486 int minidx, maxidx;
1487 unsigned sizeint[3], sizesmall[3], bitsizeint[3], size3, *luip;
1488 int k, *buf1, *buf2;
1489 int smallnum, smaller, larger, i, j, is_small, is_smaller, run, prevrun;
1490 double *lfp;
1491 float float_prec, lf,tmpdata[30];
1492 int tmp, tmpsum, *thiscoord, prevcoord[3];
1493 unsigned int tmpcoord[30];
1494 int errval=1;
1495 unsigned int bitsize;
1496
1497 bitsizeint[0] = 0;
1498 bitsizeint[1] = 0;
1499 bitsizeint[2] = 0;
1500
1501 if(xfp==NULL)
1502 return -1;
1503 size3=3*size;
1504 if(size3>xfp->buf1size) {
1505 if((xfp->buf1=(int *)malloc(sizeof(int)*size3))==NULL) {
1506 fprintf(stderr,"Cannot allocate memory for compressing coordinates.\n");
1507 return -1;
1508 }
1509 xfp->buf1size=size3;
1510 xfp->buf2size=size3*1.2;
1511 if((xfp->buf2=(int *)malloc(sizeof(int)*xfp->buf2size))==NULL) {
1512 fprintf(stderr,"Cannot allocate memory for compressing coordinates.\n");
1513 return -1;
1514 }
1515 }
1516 if(xdrfile_write_int(&size,1,xfp)==0)
1517 return -1; /* return if we could not write size */
1518 /* Dont bother with compression for three atoms or less */
1519 if(size<=9) {
1520 for(i=0;i<9*3;i++)
1521 tmpdata[i]=ptr[i];
1522 return xdrfile_write_float(tmpdata,size3,xfp)/3;
1523 /* return number of coords, not floats */
1524 }
1525 /* Compression-time if we got here. Write precision first */
1526 if (precision <= 0)
1527 precision = 1000;
1528 float_prec=precision;
1529 xdrfile_write_float(&float_prec,1,xfp);
1530 /* avoid repeated pointer dereferencing. */
1531 buf1=xfp->buf1;
1532 buf2=xfp->buf2;
1533 /* buf2[0-2] are special and do not contain actual data */
1534 buf2[0] = buf2[1] = buf2[2] = 0;
1535 minint[0] = minint[1] = minint[2] = INT_MAX;
1536 maxint[0] = maxint[1] = maxint[2] = INT_MIN;
1537 prevrun = -1;
1538 lfp = ptr;
1539 lip = buf1;
1540 mindiff = INT_MAX;
1541 oldlint1 = oldlint2 = oldlint3 = 0;
1542 while(lfp < ptr + size3 ) {
1543 /* find nearest integer */
1544 if (*lfp >= 0.0)
1545 lf = (float)*lfp * float_prec + 0.5;
1546 else
1547 lf = (float)*lfp * float_prec - 0.5;
1548 if (fabs(lf) > INT_MAX-2) {
1549 /* scaling would cause overflow */
1550 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1551 errval=0;
1552 }
1553 lint1 = lf;
1554 if (lint1 < minint[0]) minint[0] = lint1;
1555 if (lint1 > maxint[0]) maxint[0] = lint1;
1556 *lip++ = lint1;
1557 lfp++;
1558 if (*lfp >= 0.0)
1559 lf = (float)*lfp * float_prec + 0.5;
1560 else
1561 lf = (float)*lfp * float_prec - 0.5;
1562 if (fabs(lf) > INT_MAX-2) {
1563 /* scaling would cause overflow */
1564 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1565 errval=0;
1566 }
1567 lint2 = lf;
1568 if (lint2 < minint[1]) minint[1] = lint2;
1569 if (lint2 > maxint[1]) maxint[1] = lint2;
1570 *lip++ = lint2;
1571 lfp++;
1572 if (*lfp >= 0.0)
1573 lf = (float)*lfp * float_prec + 0.5;
1574 else
1575 lf = (float)*lfp * float_prec - 0.5;
1576 if (fabs(lf) > INT_MAX-2) {
1577 errval=0;
1578 }
1579 lint3 = lf;
1580 if (lint3 < minint[2]) minint[2] = lint3;
1581 if (lint3 > maxint[2]) maxint[2] = lint3;
1582 *lip++ = lint3;
1583 lfp++;
1584 diff = abs(oldlint1-lint1)+abs(oldlint2-lint2)+abs(oldlint3-lint3);
1585 if (diff < mindiff && lfp > ptr + 3)
1586 mindiff = diff;
1587 oldlint1 = lint1;
1588 oldlint2 = lint2;
1589 oldlint3 = lint3;
1590 }
1591 xdrfile_write_int(minint,3,xfp);
1592 xdrfile_write_int(maxint,3,xfp);
1593
1594 if ((float)maxint[0] - (float)minint[0] >= INT_MAX-2 ||
1595 (float)maxint[1] - (float)minint[1] >= INT_MAX-2 ||
1596 (float)maxint[2] - (float)minint[2] >= INT_MAX-2) {
1597 /* turning value in unsigned by subtracting minint
1598 * would cause overflow
1599 */
1600 fprintf(stderr,"Internal overflow compressing coordinates.\n");
1601 errval=0;
1602 }
1603 sizeint[0] = maxint[0] - minint[0]+1;
1604 sizeint[1] = maxint[1] - minint[1]+1;
1605 sizeint[2] = maxint[2] - minint[2]+1;
1606
1607 /* check if one of the sizes is to big to be multiplied */
1608 if ((sizeint[0] | sizeint[1] | sizeint[2] ) > 0xffffff) {
1609 bitsizeint[0] = sizeofint(sizeint[0]);
1610 bitsizeint[1] = sizeofint(sizeint[1]);
1611 bitsizeint[2] = sizeofint(sizeint[2]);
1612 bitsize = 0; /* flag the use of large sizes */
1613 } else {
1614 bitsize = sizeofints(3, sizeint);
1615 }
1616 lip = buf1;
1617 luip = (unsigned int *) buf1;
1618 smallidx = FIRSTIDX;
1619 while (smallidx < LASTIDX && magicints[smallidx] < mindiff) {
1620 smallidx++;
1621 }
1622 xdrfile_write_int(&smallidx,1,xfp);
1623 tmp=smallidx+8;
1624 maxidx = (LASTIDX<tmp) ? LASTIDX : tmp;
1625 minidx = maxidx - 8; /* often this equal smallidx */
1626 tmp=smallidx-1;
1627 tmp= (FIRSTIDX>tmp) ? FIRSTIDX : tmp;
1628 smaller = magicints[tmp] / 2;
1629 smallnum = magicints[smallidx] / 2;
1630 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx];
1631 larger = magicints[maxidx] / 2;
1632 i = 0;
1633 while (i < size) {
1634 is_small = 0;
1635 thiscoord = (int *)(luip) + i * 3;
1636 if (smallidx < maxidx && i >= 1 &&
1637 abs(thiscoord[0] - prevcoord[0]) < larger &&
1638 abs(thiscoord[1] - prevcoord[1]) < larger &&
1639 abs(thiscoord[2] - prevcoord[2]) < larger) {
1640 is_smaller = 1;
1641 } else if (smallidx > minidx) {
1642 is_smaller = -1;
1643 } else {
1644 is_smaller = 0;
1645 }
1646 if (i + 1 < size) {
1647 if (abs(thiscoord[0] - thiscoord[3]) < smallnum &&
1648 abs(thiscoord[1] - thiscoord[4]) < smallnum &&
1649 abs(thiscoord[2] - thiscoord[5]) < smallnum) {
1650 /* interchange first with second atom for better
1651 * compression of water molecules
1652 */
1653 tmp = thiscoord[0]; thiscoord[0] = thiscoord[3];
1654 thiscoord[3] = tmp;
1655 tmp = thiscoord[1]; thiscoord[1] = thiscoord[4];
1656 thiscoord[4] = tmp;
1657 tmp = thiscoord[2]; thiscoord[2] = thiscoord[5];
1658 thiscoord[5] = tmp;
1659 is_small = 1;
1660 }
1661 }
1662 tmpcoord[0] = thiscoord[0] - minint[0];
1663 tmpcoord[1] = thiscoord[1] - minint[1];
1664 tmpcoord[2] = thiscoord[2] - minint[2];
1665 if (bitsize == 0) {
1666 encodebits(buf2, bitsizeint[0], tmpcoord[0]);
1667 encodebits(buf2, bitsizeint[1], tmpcoord[1]);
1668 encodebits(buf2, bitsizeint[2], tmpcoord[2]);
1669 } else {
1670 encodeints(buf2, 3, bitsize, sizeint, tmpcoord);
1671 }
1672 prevcoord[0] = thiscoord[0];
1673 prevcoord[1] = thiscoord[1];
1674 prevcoord[2] = thiscoord[2];
1675 thiscoord = thiscoord + 3;
1676 i++;
1677
1678 run = 0;
1679 if (is_small == 0 && is_smaller == -1)
1680 is_smaller = 0;
1681 while (is_small && run < 8*3) {
1682 tmpsum=0;
1683 for(j=0;j<3;j++) {
1684 tmp=thiscoord[j] - prevcoord[j];
1685 tmpsum+=tmp*tmp;
1686 }
1687 if (is_smaller == -1 && tmpsum >= smaller * smaller) {
1688 is_smaller = 0;
1689 }
1690
1691 tmpcoord[run++] = thiscoord[0] - prevcoord[0] + smallnum;
1692 tmpcoord[run++] = thiscoord[1] - prevcoord[1] + smallnum;
1693 tmpcoord[run++] = thiscoord[2] - prevcoord[2] + smallnum;
1694
1695 prevcoord[0] = thiscoord[0];
1696 prevcoord[1] = thiscoord[1];
1697 prevcoord[2] = thiscoord[2];
1698
1699 i++;
1700 thiscoord = thiscoord + 3;
1701 is_small = 0;
1702 if (i < size &&
1703 abs(thiscoord[0] - prevcoord[0]) < smallnum &&
1704 abs(thiscoord[1] - prevcoord[1]) < smallnum &&
1705 abs(thiscoord[2] - prevcoord[2]) < smallnum) {
1706 is_small = 1;
1707 }
1708 }
1709 if (run != prevrun || is_smaller != 0) {
1710 prevrun = run;
1711 encodebits(buf2, 1, 1); /* flag the change in run-length */
1712 encodebits(buf2, 5, run+is_smaller+1);
1713 } else {
1714 encodebits(buf2, 1, 0); /* flag the fact that runlength did not change */
1715 }
1716 for (k=0; k < run; k+=3) {
1717 encodeints(buf2, 3, smallidx, sizesmall, &tmpcoord[k]);
1718 }
1719 if (is_smaller != 0) {
1720 smallidx += is_smaller;
1721 if (is_smaller < 0) {
1722 smallnum = smaller;
1723 smaller = magicints[smallidx-1] / 2;
1724 } else {
1725 smaller = smallnum;
1726 smallnum = magicints[smallidx] / 2;
1727 }
1728 sizesmall[0] = sizesmall[1] = sizesmall[2] = magicints[smallidx];
1729 }
1730 }
1731 if (buf2[1] != 0) buf2[0]++;
1732 xdrfile_write_int(buf2,1,xfp); /* buf2[0] holds the length in bytes */
1733 tmp=xdrfile_write_opaque((char *)&(buf2[3]),(unsigned int)buf2[0],xfp);
1734 if(tmp==(unsigned int)buf2[0])
1735 return size;
1736 else
1737 return -1;
1738 }
1739
1740
1741 /* Dont try do document Fortran interface, since
1742 * Doxygen barfs at the F77_FUNC macro
1743 */
1744 #ifndef DOXYGEN
1745
1746 /*************************************************************
1747 * Fortran77 interface for reading/writing portable data *
1748 * The routine are not threadsafe when called from Fortran *
1749 * (as they are when called from C) unless you compile with *
1750 * this file with posix thread support. *
1751 * Note that these are not multithread-safe. *
1752 *************************************************************/
1753 #define MAX_FORTRAN_XDR 1024
1754 static XDRFILE *f77xdr[MAX_FORTRAN_XDR]; /* array of file handles */
1755 static int f77init = 1; /* zero array first time */
1756
1757 /* internal to this file: C<-->Fortran string conversion */
1758 static int ftocstr(char *dest, int dest_len, char *src, int src_len);
1759 static int ctofstr(char *dest, int dest_len, char *src);
1760
1761
1762 void
F77_FUNC(xdropen,XDROPEN)1763 F77_FUNC(xdropen,XDROPEN)(int *fid, char *filename, char *mode,
1764 int fn_len, int mode_len)
1765 {
1766 char cfilename[512];
1767 char cmode[5];
1768 int i;
1769
1770 /* zero array at first invocation */
1771 if(f77init) {
1772 for(i=0;i<MAX_FORTRAN_XDR;i++)
1773 f77xdr[i]=NULL;
1774 f77init=0;
1775 }
1776 i=0;
1777
1778 /* nf77xdr is always smaller or equal to MAX_FORTRAN_XDR */
1779 while(i<MAX_FORTRAN_XDR && f77xdr[i]!=NULL)
1780 i++;
1781 if(i==MAX_FORTRAN_XDR) {
1782 *fid = -1;
1783 } else if (ftocstr(cfilename, sizeof(cfilename), filename, fn_len)) {
1784 *fid = -1;
1785 } else if (ftocstr(cmode, sizeof(cmode), mode,mode_len)) {
1786 *fid = -1;
1787 } else {
1788 f77xdr[i]=xdrfile_open(cfilename,cmode);
1789 /* return the index in the array as a fortran file handle */
1790 *fid=i;
1791 }
1792 }
1793
1794 void
F77_FUNC(xdrclose,XDRCLOSE)1795 F77_FUNC(xdrclose,XDRCLOSE)(int *fid)
1796 {
1797 /* first close it */
1798 xdrfile_close(f77xdr[*fid]);
1799 /* the remove it from file handle list */
1800 f77xdr[*fid]=NULL;
1801 }
1802
1803
1804 void
F77_FUNC(xdrrint,XDRRINT)1805 F77_FUNC(xdrrint,XDRRINT)(int *fid, int *data, int *ndata, int *ret)
1806 {
1807 *ret = xdrfile_read_int(data,*ndata,f77xdr[*fid]);
1808 }
1809
1810 void
F77_FUNC(xdrwint,XDRWINT)1811 F77_FUNC(xdrwint,XDRWINT)(int *fid, int *data, int *ndata, int *ret)
1812 {
1813 *ret = xdrfile_write_int(data,*ndata,f77xdr[*fid]);
1814 }
1815
1816 void
F77_FUNC(xdrruint,XDRRUINT)1817 F77_FUNC(xdrruint,XDRRUINT)(int *fid, unsigned int *data, int *ndata, int *ret)
1818 {
1819 *ret = xdrfile_read_uint(data,*ndata,f77xdr[*fid]);
1820 }
1821
1822 void
F77_FUNC(xdrwuint,XDRWUINT)1823 F77_FUNC(xdrwuint,XDRWUINT)(int *fid, unsigned int *data, int *ndata, int *ret)
1824 {
1825 *ret = xdrfile_write_uint(data,*ndata,f77xdr[*fid]);
1826 }
1827
1828 void
F77_FUNC(xdrrchar,XDRRCHAR)1829 F77_FUNC(xdrrchar,XDRRCHAR)(int *fid, char *ip, int *ndata, int *ret)
1830 {
1831 *ret = xdrfile_read_char(ip,*ndata,f77xdr[*fid]);
1832 }
1833
1834 void
F77_FUNC(xdrwchar,XDRWCHAR)1835 F77_FUNC(xdrwchar,XDRWCHAR)(int *fid, char *ip, int *ndata, int *ret)
1836 {
1837 *ret = xdrfile_write_char(ip,*ndata,f77xdr[*fid]);
1838 }
1839
1840 void
F77_FUNC(xdrruchar,XDRRUCHAR)1841 F77_FUNC(xdrruchar,XDRRUCHAR)(int *fid, unsigned char *ip, int *ndata, int *ret)
1842 {
1843 *ret = xdrfile_read_uchar(ip,*ndata,f77xdr[*fid]);
1844 }
1845
1846 void
F77_FUNC(xdrwuchar,XDRWUCHAR)1847 F77_FUNC(xdrwuchar,XDRWUCHAR)(int *fid, unsigned char *ip, int *ndata, int *ret)
1848 {
1849 *ret = xdrfile_write_uchar(ip,*ndata,f77xdr[*fid]);
1850 }
1851
1852 void
F77_FUNC(xdrrshort,XDRRSHORT)1853 F77_FUNC(xdrrshort,XDRRSHORT)(int *fid, short *ip, int *ndata, int *ret)
1854 {
1855 *ret = xdrfile_read_short(ip,*ndata,f77xdr[*fid]);
1856 }
1857
1858 void
F77_FUNC(xdrwshort,XDRWSHORT)1859 F77_FUNC(xdrwshort,XDRWSHORT)(int *fid, short *ip, int *ndata, int *ret)
1860 {
1861 *ret = xdrfile_write_short(ip,*ndata,f77xdr[*fid]);
1862 }
1863
1864 void
F77_FUNC(xdrrushort,XDRRUSHORT)1865 F77_FUNC(xdrrushort,XDRRUSHORT)(int *fid, unsigned short *ip, int *ndata, int *ret)
1866 {
1867 *ret = xdrfile_read_ushort(ip,*ndata,f77xdr[*fid]);
1868 }
1869
1870 void
F77_FUNC(xdrwushort,XDRWUSHORT)1871 F77_FUNC(xdrwushort,XDRWUSHORT)(int *fid, unsigned short *ip, int *ndata, int *ret)
1872 {
1873 *ret = xdrfile_write_ushort(ip,*ndata,f77xdr[*fid]);
1874 }
1875
1876 void
F77_FUNC(xdrrsingle,XDRRSINGLE)1877 F77_FUNC(xdrrsingle,XDRRSINGLE)(int *fid, float *data, int *ndata, int *ret)
1878 {
1879 *ret = xdrfile_read_float(data,*ndata,f77xdr[*fid]);
1880 }
1881
1882 void
F77_FUNC(xdrwsingle,XDRWSINGLE)1883 F77_FUNC(xdrwsingle,XDRWSINGLE)(int *fid, float *data, int *ndata, int *ret)
1884 {
1885 *ret = xdrfile_write_float(data,*ndata,f77xdr[*fid]);
1886 }
1887
1888 void
F77_FUNC(xdrrdouble,XDRRDOUBLE)1889 F77_FUNC(xdrrdouble,XDRRDOUBLE)(int *fid, double *data, int *ndata, int *ret)
1890 {
1891 *ret = xdrfile_read_double(data,*ndata,f77xdr[*fid]);
1892 }
1893
1894 void
F77_FUNC(xdrwdouble,XDRWDOUBLE)1895 F77_FUNC(xdrwdouble,XDRWDOUBLE)(int *fid, double *data, int *ndata, int *ret)
1896 {
1897 *ret = xdrfile_write_double(data,*ndata,f77xdr[*fid]);
1898 }
1899
ftocstr(char * dest,int destlen,char * src,int srclen)1900 static int ftocstr(char *dest, int destlen, char *src, int srclen)
1901 {
1902 char *p;
1903
1904 p = src + srclen;
1905 while ( --p >= src && *p == ' ' );
1906 srclen = p - src + 1;
1907 destlen--;
1908 dest[0] = 0;
1909 if (srclen > destlen)
1910 return 1;
1911 while (srclen--)
1912 (*dest++ = *src++);
1913 *dest = '\0';
1914 return 0;
1915 }
1916
1917
ctofstr(char * dest,int destlen,char * src)1918 static int ctofstr(char *dest, int destlen, char *src)
1919 {
1920 while (destlen && *src) {
1921 *dest++ = *src++;
1922 destlen--;
1923 }
1924 while (destlen--)
1925 *dest++ = ' ';
1926 return 0;
1927 }
1928
1929
1930 void
F77_FUNC(xdrrstring,XDRRSTRING)1931 F77_FUNC(xdrrstring,XDRRSTRING)(int *fid, char *str, int *ret, int len)
1932 {
1933 char *cstr;
1934
1935 if((cstr=(char*)malloc((len+1)*sizeof(char)))==NULL) {
1936 *ret = 0;
1937 return;
1938 }
1939 if (ftocstr(cstr, len+1, str, len)) {
1940 *ret = 0;
1941 free(cstr);
1942 return;
1943 }
1944
1945 *ret = xdrfile_read_string(cstr, len+1,f77xdr[*fid]);
1946 ctofstr( str, len , cstr);
1947 free(cstr);
1948 }
1949
1950 void
F77_FUNC(xdrwstring,XDRWSTRING)1951 F77_FUNC(xdrwstring,XDRWSTRING)(int *fid, char *str, int *ret, int len)
1952 {
1953 char *cstr;
1954
1955 if((cstr=(char*)malloc((len+1)*sizeof(char)))==NULL) {
1956 *ret = 0;
1957 return;
1958 }
1959 if (ftocstr(cstr, len+1, str, len)) {
1960 *ret = 0;
1961 free(cstr);
1962 return;
1963 }
1964
1965 *ret = xdrfile_write_string(cstr, f77xdr[*fid]);
1966 ctofstr( str, len , cstr);
1967 free(cstr);
1968 }
1969
1970 void
F77_FUNC(xdrropaque,XDRROPAQUE)1971 F77_FUNC(xdrropaque,XDRROPAQUE)(int *fid, char *data, int *ndata, int *ret)
1972 {
1973 *ret = xdrfile_read_opaque(data,*ndata,f77xdr[*fid]);
1974 }
1975
1976 void
F77_FUNC(xdrwopaque,XDRWOPAQUE)1977 F77_FUNC(xdrwopaque,XDRWOPAQUE)(int *fid, char *data, int *ndata, int *ret)
1978 {
1979 *ret = xdrfile_write_opaque(data,*ndata,f77xdr[*fid]);
1980 }
1981
1982
1983 /* Write single-precision compressed 3d coordinates */
1984 void
F77_FUNC(xdrccs,XDRCCS)1985 F77_FUNC(xdrccs,XDRCCS)(int *fid, float *data, int *ncoord,
1986 float *precision, int *ret)
1987 {
1988 *ret = xdrfile_compress_coord_float(data,*ncoord,*precision,f77xdr[*fid]);
1989 }
1990
1991
1992 /* Read single-precision compressed 3d coordinates */
1993 void
F77_FUNC(xdrdcs,XDRDCS)1994 F77_FUNC(xdrdcs,XDRDCS)(int *fid, float *data, int *ncoord,
1995 float *precision, int *ret)
1996 {
1997 *ret = xdrfile_decompress_coord_float(data,ncoord,precision,f77xdr[*fid]);
1998 }
1999
2000
2001 /* Write compressed 3d coordinates from double precision data */
2002 void
F77_FUNC(xdrccd,XDRCCD)2003 F77_FUNC(xdrccd,XDRCCD)(int *fid, double *data, int *ncoord,
2004 double *precision, int *ret)
2005 {
2006 *ret = xdrfile_compress_coord_double(data,*ncoord,*precision,f77xdr[*fid]);
2007 }
2008
2009 /* Read compressed 3d coordinates into double precision data */
2010 void
F77_FUNC(xddcd,XDRDCD)2011 F77_FUNC(xddcd,XDRDCD)(int *fid, double *data, int *ncoord,
2012 double *precision, int *ret)
2013 {
2014 *ret = xdrfile_decompress_coord_double(data,ncoord,precision,f77xdr[*fid]);
2015 }
2016
2017
2018
2019
2020
2021
2022
2023 #endif /* DOXYGEN */
2024
2025 /*************************************************************
2026 * End of higher-level routines - dont change things below! *
2027 *************************************************************/
2028
2029
2030
2031
2032
2033
2034
2035
2036
2037
2038
2039
2040
2041
2042
2043
2044
2045
2046
2047
2048 /*************************************************************
2049 * The rest of this file contains our own implementation *
2050 * of the XDR calls in case you are compiling without them. *
2051 * You do NOT want to change things here since it would make *
2052 * things incompatible with the standard RPC/XDR routines. *
2053 *************************************************************/
2054 #ifndef HAVE_RPC_XDR_H
2055
2056 /*
2057 * What follows is a modified version of the Sun XDR code. For reference
2058 * we include their copyright and license:
2059 *
2060 * Sun RPC is a product of Sun Microsystems, Inc. and is provided for
2061 * unrestricted use provided that this legend is included on all tape
2062 * media and as a part of the software program in whole or part. Users
2063 * may copy or modify Sun RPC without charge, but are not authorized
2064 * to license or distribute it to anyone else except as part of a product or
2065 * program developed by the user.
2066 *
2067 * SUN RPC IS PROVIDED AS IS WITH NO WARRANTIES OF ANY KIND INCLUDING THE
2068 * WARRANTIES OF DESIGN, MERCHANTIBILITY AND FITNESS FOR A PARTICULAR
2069 * PURPOSE, OR ARISING FROM A COURSE OF DEALING, USAGE OR TRADE PRACTICE.
2070 *
2071 * Sun RPC is provided with no support and without any obligation on the
2072 * part of Sun Microsystems, Inc. to assist in its use, correction,
2073 * modification or enhancement.
2074 *
2075 * SUN MICROSYSTEMS, INC. SHALL HAVE NO LIABILITY WITH RESPECT TO THE
2076 * INFRINGEMENT OF COPYRIGHTS, TRADE SECRETS OR ANY PATENTS BY SUN RPC
2077 * OR ANY PART THEREOF.
2078 *
2079 * In no event will Sun Microsystems, Inc. be liable for any lost revenue
2080 * or profits or other special, indirect and consequential damages, even if
2081 * Sun has been advised of the possibility of such damages.
2082 *
2083 * Sun Microsystems, Inc.
2084 * 2550 Garcia Avenue
2085 * Mountain View, California 94043
2086 */
2087
2088 /* INT_MAX is defined in limits.h according to ANSI C */
2089 #if (INT_MAX > 2147483647)
2090 # error Error: Cannot use builtin XDR support when size of int
2091 # error is larger than 4 bytes. Use your system XDR libraries
2092 # error instead, or modify the source code in xdrfile.c
2093 #endif /* Check for 4 byte int type */
2094
2095
2096
2097
2098
2099 typedef int (*xdrproc_t) (XDR *, void *,...);
2100
2101 #define xdr_getlong(xdrs, longp) \
2102 (*(xdrs)->x_ops->x_getlong)(xdrs, longp)
2103 #define xdr_putlong(xdrs, longp) \
2104 (*(xdrs)->x_ops->x_putlong)(xdrs, longp)
2105 #define xdr_getbytes(xdrs, addr, len) \
2106 (*(xdrs)->x_ops->x_getbytes)(xdrs, addr, len)
2107 #define xdr_putbytes(xdrs, addr, len) \
2108 (*(xdrs)->x_ops->x_putbytes)(xdrs, addr, len)
2109
2110 #define BYTES_PER_XDR_UNIT 4
2111 static char xdr_zero[BYTES_PER_XDR_UNIT] = {0, 0, 0, 0};
2112
2113 static int32_t
xdr_swapbytes(int32_t x)2114 xdr_swapbytes(int32_t x)
2115 {
2116 int32_t y,i;
2117 char *px=(char *)&x;
2118 char *py=(char *)&y;
2119
2120 for(i=0;i<4;i++)
2121 py[i]=px[3-i];
2122
2123 return y;
2124 }
2125
2126 static int32_t
xdr_htonl(int32_t x)2127 xdr_htonl(int32_t x)
2128 {
2129 int s=0x1234;
2130 if( *((char *)&s)==(char)0x34) {
2131 /* smallendian,swap bytes */
2132 return xdr_swapbytes(x);
2133 } else {
2134 /* bigendian, do nothing */
2135 return x;
2136 }
2137 }
2138
2139 static int32_t
xdr_ntohl(int x)2140 xdr_ntohl(int x)
2141 {
2142 int s=0x1234;
2143 if( *((char *)&s)==(char)0x34) {
2144 /* smallendian, swap bytes */
2145 return xdr_swapbytes(x);
2146 } else {
2147 /* bigendian, do nothing */
2148 return x;
2149 }
2150 }
2151
2152 static int
xdr_int(XDR * xdrs,int * ip)2153 xdr_int (XDR *xdrs, int *ip)
2154 {
2155 int32_t i32;
2156
2157 switch (xdrs->x_op)
2158 {
2159 case XDR_ENCODE:
2160 i32 = (int32_t) *ip;
2161 return xdr_putlong (xdrs, &i32);
2162
2163 case XDR_DECODE:
2164 if (!xdr_getlong (xdrs, &i32))
2165 {
2166 return 0;
2167 }
2168 *ip = (int) i32;
2169 case XDR_FREE:
2170 return 1;
2171 }
2172 return 0;
2173 }
2174
2175 static int
xdr_u_int(XDR * xdrs,unsigned int * up)2176 xdr_u_int (XDR *xdrs, unsigned int *up)
2177 {
2178 uint32_t ui32;
2179
2180 switch (xdrs->x_op)
2181 {
2182 case XDR_ENCODE:
2183 ui32 = (uint32_t) * up;
2184 return xdr_putlong (xdrs, (int32_t *)&ui32);
2185
2186 case XDR_DECODE:
2187 if (!xdr_getlong (xdrs, (int32_t *)&ui32))
2188 {
2189 return 0;
2190 }
2191 *up = (uint32_t) ui32;
2192 case XDR_FREE:
2193 return 1;
2194 }
2195 return 0;
2196 }
2197
2198 static int
xdr_short(XDR * xdrs,short * sp)2199 xdr_short (XDR *xdrs, short *sp)
2200 {
2201 int32_t i32;
2202
2203 switch (xdrs->x_op)
2204 {
2205 case XDR_ENCODE:
2206 i32 = (int32_t) *sp;
2207 return xdr_putlong (xdrs, &i32);
2208
2209 case XDR_DECODE:
2210 if (!xdr_getlong (xdrs, &i32))
2211 {
2212 return 0;
2213 }
2214 *sp = (short) i32;
2215 return 1;
2216
2217 case XDR_FREE:
2218 return 1;
2219 }
2220 return 0;
2221 }
2222
2223 static int
xdr_u_short(XDR * xdrs,unsigned short * sp)2224 xdr_u_short (XDR *xdrs, unsigned short *sp)
2225 {
2226 uint32_t ui32;
2227
2228 switch (xdrs->x_op)
2229 {
2230 case XDR_ENCODE:
2231 ui32 = (uint32_t) *sp;
2232 return xdr_putlong (xdrs, (int32_t *)&ui32);
2233
2234 case XDR_DECODE:
2235 if (!xdr_getlong (xdrs, (int32_t *)&ui32))
2236 {
2237 return 0;
2238 }
2239 *sp = (unsigned short) ui32;
2240 return 1;
2241
2242 case XDR_FREE:
2243 return 1;
2244 }
2245 return 0;
2246 }
2247
2248 static int
xdr_char(XDR * xdrs,char * cp)2249 xdr_char (XDR *xdrs, char *cp)
2250 {
2251 int i;
2252
2253 i = (*cp);
2254 if (!xdr_int (xdrs, &i))
2255 {
2256 return 0;
2257 }
2258 *cp = i;
2259 return 1;
2260 }
2261
2262 static int
xdr_u_char(XDR * xdrs,unsigned char * cp)2263 xdr_u_char (XDR *xdrs, unsigned char *cp)
2264 {
2265 unsigned int u;
2266
2267 u = (*cp);
2268 if (!xdr_u_int (xdrs, &u))
2269 {
2270 return 0;
2271 }
2272 *cp = u;
2273 return 1;
2274 }
2275
2276 /*
2277 * XDR opaque data
2278 * Allows the specification of a fixed size sequence of opaque bytes.
2279 * cp points to the opaque object and cnt gives the byte length.
2280 */
2281 static int
xdr_opaque(XDR * xdrs,char * cp,unsigned int cnt)2282 xdr_opaque (XDR *xdrs, char *cp, unsigned int cnt)
2283 {
2284 unsigned int rndup;
2285 static char crud[BYTES_PER_XDR_UNIT];
2286
2287 /*
2288 * if no data we are done
2289 */
2290 if (cnt == 0)
2291 return 1;
2292
2293 /*
2294 * round byte count to full xdr units
2295 */
2296 rndup = cnt % BYTES_PER_XDR_UNIT;
2297 if (rndup > 0)
2298 rndup = BYTES_PER_XDR_UNIT - rndup;
2299
2300 switch (xdrs->x_op)
2301 {
2302 case XDR_DECODE:
2303 if (!xdr_getbytes (xdrs, cp, cnt))
2304 {
2305 return 0;
2306 }
2307 if (rndup == 0)
2308 return 1;
2309 return xdr_getbytes (xdrs, (char *)crud, rndup);
2310
2311 case XDR_ENCODE:
2312 if (!xdr_putbytes (xdrs, cp, cnt))
2313 {
2314 return 0;
2315 }
2316 if (rndup == 0)
2317 return 1;
2318 return xdr_putbytes (xdrs, xdr_zero, rndup);
2319
2320 case XDR_FREE:
2321 return 1;
2322 }
2323 #undef BYTES_PER_XDR_UNIT
2324 return 0;
2325 }
2326
2327
2328 /*
2329 * XDR null terminated ASCII strings
2330 */
2331 static int
xdr_string(XDR * xdrs,char ** cpp,unsigned int maxsize)2332 xdr_string (XDR *xdrs, char **cpp, unsigned int maxsize)
2333 {
2334 char *sp = *cpp; /* sp is the actual string pointer */
2335 unsigned int size;
2336 unsigned int nodesize;
2337
2338 /*
2339 * first deal with the length since xdr strings are counted-strings
2340 */
2341 switch (xdrs->x_op)
2342 {
2343 case XDR_FREE:
2344 if (sp == NULL)
2345 {
2346 return 1; /* already free */
2347 }
2348 /* fall through... */
2349 case XDR_ENCODE:
2350 if (sp == NULL)
2351 return 0;
2352 size = strlen (sp);
2353 break;
2354 case XDR_DECODE:
2355 break;
2356 }
2357 if (!xdr_u_int (xdrs, &size))
2358 {
2359 return 0;
2360 }
2361 if (size > maxsize)
2362 {
2363 return 0;
2364 }
2365 nodesize = size + 1;
2366
2367 /*
2368 * now deal with the actual bytes
2369 */
2370 switch (xdrs->x_op)
2371 {
2372 case XDR_DECODE:
2373 if (nodesize == 0)
2374 {
2375 return 1;
2376 }
2377 if (sp == NULL)
2378 *cpp = sp = (char *) malloc (nodesize);
2379 if (sp == NULL)
2380 {
2381 (void) fputs ("xdr_string: out of memory\n", stderr);
2382 return 0;
2383 }
2384 sp[size] = 0;
2385 /* fall into ... */
2386
2387 case XDR_ENCODE:
2388 return xdr_opaque (xdrs, sp, size);
2389
2390 case XDR_FREE:
2391 free (sp);
2392 *cpp = NULL;
2393 return 1;
2394 }
2395 return 0;
2396 }
2397
2398
2399
2400 /* Floating-point stuff */
2401
2402 static int
xdr_float(XDR * xdrs,float * fp)2403 xdr_float(XDR *xdrs, float *fp)
2404 {
2405 switch (xdrs->x_op) {
2406
2407 case XDR_ENCODE:
2408 if (sizeof(float) == sizeof(int32_t))
2409 return (xdr_putlong(xdrs, (int32_t *)fp));
2410 else if (sizeof(float) == sizeof(int)) {
2411 int32_t tmp = *(int *)fp;
2412 return (xdr_putlong(xdrs, &tmp));
2413 }
2414 break;
2415
2416 case XDR_DECODE:
2417 if (sizeof(float) == sizeof(int32_t))
2418 return (xdr_getlong(xdrs, (int32_t *)fp));
2419 else if (sizeof(float) == sizeof(int)) {
2420 int32_t tmp;
2421 if (xdr_getlong(xdrs, &tmp)) {
2422 *(int *)fp = tmp;
2423 return (1);
2424 }
2425 }
2426 break;
2427
2428 case XDR_FREE:
2429 return (1);
2430 }
2431 return (0);
2432 }
2433
2434
2435 static int
xdr_double(XDR * xdrs,double * dp)2436 xdr_double(XDR *xdrs, double *dp)
2437 {
2438 /* Gromacs detects floating-point stuff at compile time, which is faster */
2439 #ifdef GROMACS
2440 # ifndef FLOAT_FORMAT_IEEE754
2441 # error non-IEEE floating point system, or you defined GROMACS yourself...
2442 # endif
2443 int LSW;
2444 # ifdef IEEE754_BIG_ENDIAN_WORD_ORDER
2445 int LSW=1;
2446 # else
2447 int LSW=0;
2448 # endif /* Big endian word order */
2449 #else
2450 /* Outside Gromacs we rely on dynamic detection of FP order. */
2451 int LSW; /* Least significant fp word */
2452
2453 double x=0.987654321; /* Just a number */
2454 unsigned char ix = *((char *)&x);
2455
2456 /* Possible representations in IEEE double precision:
2457 * (S=small endian, B=big endian)
2458 *
2459 * Byte order, Word order, Hex
2460 * S S b8 56 0e 3c dd 9a ef 3f
2461 * B S 3c 0e 56 b8 3f ef 9a dd
2462 * S B dd 9a ef 3f b8 56 0e 3c
2463 * B B 3f ef 9a dd 3c 0e 56 b8
2464 */
2465 if(ix==0xdd || ix==0x3f)
2466 LSW=1; /* Big endian word order */
2467 else if(ix==0xb8 || ix==0x3c)
2468 LSW=0; /* Small endian word order */
2469 else { /* Catch strange errors */
2470 fprintf(stderr,"Cannot detect floating-point word order.\n"
2471 "Do you have a non-IEEE system?\n"
2472 "Use system XDR libraries or fix xdr_double().\n");
2473 abort();
2474 }
2475 #endif /* end of dynamic detection of fp word order */
2476
2477 switch (xdrs->x_op) {
2478
2479 case XDR_ENCODE:
2480 if (2*sizeof(int32_t) == sizeof(double)) {
2481 int32_t *lp = (int32_t *)dp;
2482 return (xdr_putlong(xdrs, lp+!LSW) &&
2483 xdr_putlong(xdrs, lp+LSW));
2484 } else if (2*sizeof(int) == sizeof(double)) {
2485 int *ip = (int *)dp;
2486 int32_t tmp[2];
2487 tmp[0] = ip[!LSW];
2488 tmp[1] = ip[LSW];
2489 return (xdr_putlong(xdrs, tmp) &&
2490 xdr_putlong(xdrs, tmp+1));
2491 }
2492 break;
2493
2494 case XDR_DECODE:
2495 if (2*sizeof(int32_t) == sizeof(double)) {
2496 int32_t *lp = (int32_t *)dp;
2497 return (xdr_getlong(xdrs, lp+!LSW) &&
2498 xdr_getlong(xdrs, lp+LSW));
2499 } else if (2*sizeof(int) == sizeof(double)) {
2500 int *ip = (int *)dp;
2501 int32_t tmp[2];
2502 if (xdr_getlong(xdrs, tmp+!LSW) &&
2503 xdr_getlong(xdrs, tmp+LSW)) {
2504 ip[0] = tmp[0];
2505 ip[1] = tmp[1];
2506 return (1);
2507 }
2508 }
2509 break;
2510
2511 case XDR_FREE:
2512 return (1);
2513 }
2514 return (0);
2515 }
2516
2517
2518 static int xdrstdio_getlong (XDR *, int32_t *);
2519 static int xdrstdio_putlong (XDR *, int32_t *);
2520 static int xdrstdio_getbytes (XDR *, char *, unsigned int);
2521 static int xdrstdio_putbytes (XDR *, char *, unsigned int);
2522 static int64_t xdrstdio_getpos (XDR *);
2523 static int xdrstdio_setpos (XDR *, int64_t, int);
2524 static void xdrstdio_destroy (XDR *);
2525
2526 /*
2527 * Ops vector for stdio type XDR
2528 */
2529 static const struct xdr_ops xdrstdio_ops =
2530 {
2531 xdrstdio_getlong, /* deserialize a long int */
2532 xdrstdio_putlong, /* serialize a long int */
2533 xdrstdio_getbytes, /* deserialize counted bytes */
2534 xdrstdio_putbytes, /* serialize counted bytes */
2535 xdrstdio_getpos, /* get offset in the stream */
2536 xdrstdio_setpos, /* set offset in the stream */
2537 xdrstdio_destroy, /* destroy stream */
2538 };
2539
2540 /*
2541 * Initialize a stdio xdr stream.
2542 * Sets the xdr stream handle xdrs for use on the stream file.
2543 * Operation flag is set to op.
2544 */
2545 static void
xdrstdio_create(XDR * xdrs,FILE * file,enum xdr_op op)2546 xdrstdio_create (XDR *xdrs, FILE *file, enum xdr_op op)
2547 {
2548 xdrs->x_op = op;
2549
2550 xdrs->x_ops = (struct xdr_ops *) &xdrstdio_ops;
2551 xdrs->x_private = (char *) file;
2552 }
2553
2554 /*
2555 * Destroy a stdio xdr stream.
2556 * Cleans up the xdr stream handle xdrs previously set up by xdrstdio_create.
2557 */
2558 static void
xdrstdio_destroy(XDR * xdrs)2559 xdrstdio_destroy (XDR *xdrs)
2560 {
2561 (void) fflush ((FILE *) xdrs->x_private);
2562 /* xx should we close the file ?? */
2563 }
2564
2565 static int
xdrstdio_getlong(XDR * xdrs,int32_t * lp)2566 xdrstdio_getlong (XDR *xdrs, int32_t *lp)
2567 {
2568 int32_t mycopy;
2569
2570 if (fread ((char *) & mycopy, 4, 1, (FILE *) xdrs->x_private) != 1)
2571 return 0;
2572 *lp = (int32_t) xdr_ntohl (mycopy);
2573 return 1;
2574 }
2575
2576 static int
xdrstdio_putlong(XDR * xdrs,int32_t * lp)2577 xdrstdio_putlong (XDR *xdrs, int32_t *lp)
2578 {
2579 int32_t mycopy = xdr_htonl (*lp);
2580 lp = &mycopy;
2581 if (fwrite ((char *) lp, 4, 1, (FILE *) xdrs->x_private) != 1)
2582 return 0;
2583 return 1;
2584 }
2585
2586 static int
xdrstdio_getbytes(XDR * xdrs,char * addr,unsigned int len)2587 xdrstdio_getbytes (XDR *xdrs, char *addr, unsigned int len)
2588 {
2589 if ((len != 0) && (fread (addr, (int) len, 1,
2590 (FILE *) xdrs->x_private) != 1))
2591 return 0;
2592 return 1;
2593 }
2594
2595 static int
xdrstdio_putbytes(XDR * xdrs,char * addr,unsigned int len)2596 xdrstdio_putbytes (XDR *xdrs, char *addr, unsigned int len)
2597 {
2598 if ((len != 0) && (fwrite (addr, (int) len, 1,
2599 (FILE *) xdrs->x_private) != 1))
2600 return 0;
2601 return 1;
2602 }
2603
2604
2605 static int64_t
xdrstdio_getpos(XDR * xdrs)2606 xdrstdio_getpos (XDR *xdrs)
2607 {
2608 #ifdef _WIN32
2609 return _ftelli64((FILE *) xdrs->x_private);
2610 #else /* __unix__ */
2611 return ftello((FILE *) xdrs->x_private);
2612 #endif
2613 }
2614
2615 static int
xdrstdio_setpos(XDR * xdrs,int64_t pos,int whence)2616 xdrstdio_setpos (XDR *xdrs, int64_t pos, int whence)
2617 {
2618 /* A reason for failure can be filesystem limits on allocation units,
2619 * before the actual off_t overflow (ext3, with a 4K clustersize,
2620 * has a 16TB limit).*/
2621 /* We return errno relying on the fact that it is never set to 0 on
2622 * success, which means that if an error occurrs it'll never be the same
2623 * as exdrOK, and xdr_seek won't be confused.*/
2624 #ifdef _WIN32
2625 return _fseeki64((FILE *) xdrs->x_private, pos, whence) < 0 ? errno : exdrOK;
2626 #else /* __unix__ */
2627 return fseeko((FILE *) xdrs->x_private, pos, whence) < 0 ? errno : exdrOK;
2628 #endif
2629
2630
2631 }
2632
2633
xdr_tell(XDRFILE * xd)2634 int64_t xdr_tell(XDRFILE *xd)
2635 /* Reads position in file */
2636 {
2637 return (int64_t)xdrstdio_getpos(xd->xdr);
2638 }
2639
xdr_seek(XDRFILE * xd,int64_t pos,int whence)2640 int xdr_seek(XDRFILE *xd, int64_t pos, int whence)
2641 /* Seeks to position in file */
2642 {
2643 int result;
2644 if ((result = xdrstdio_setpos(xd->xdr, (int64_t) pos, whence)) != 0)
2645 return result;
2646
2647 return exdrOK;
2648 }
2649
2650
2651 #endif /* HAVE_RPC_XDR_H not defined */
2652