From beaa92cb293a4147aef8ed03027500804535ed96 Mon Sep 17 00:00:00 2001 From: Magnus Lundborg Date: Mon, 21 Oct 2013 09:31:05 +0200 Subject: Fixed compiler warnings and linking errors in MSVC. Changed tabs to spaces in tng_compression functions. diff --git a/include/compression/coder.h b/include/compression/coder.h index 5cef38a..570bc6d 100644 --- a/include/compression/coder.h +++ b/include/compression/coder.h @@ -13,6 +13,14 @@ #ifndef CODER_H #define CODER_H +#ifndef DECLSPECDLLEXPORT +#ifdef USE_WINDOWS +#define DECLSPECDLLEXPORT __declspec(dllexport) +#else /* USE_WINDOWS */ +#define DECLSPECDLLEXPORT +#endif /* USE_WINDOWS */ +#endif /* DECLSPECDLLEXPORT */ + struct coder { unsigned int pack_temporary; @@ -21,22 +29,22 @@ struct coder int stat_numval; }; -struct coder *Ptngc_coder_init(void); -void Ptngc_coder_deinit(struct coder *coder); -unsigned char *Ptngc_pack_array(struct coder *coder,int *input, int *length, int coding, int coding_parameter, int natoms, int speed); -int Ptngc_unpack_array(struct coder *coder,unsigned char *packed,int *output, int length, int coding, int coding_parameter, int natoms); -unsigned char *Ptngc_pack_array_xtc2(struct coder *coder,int *input, int *length); -int Ptngc_unpack_array_xtc2(struct coder *coder,unsigned char *packed,int *output, int length); -unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int speed); -int Ptngc_unpack_array_xtc3(unsigned char *packed,int *output, int length, int natoms); - -void Ptngc_out8bits(struct coder *coder, unsigned char **output); -void Ptngc_pack_flush(struct coder *coder,unsigned char **output); -void Ptngc_write_pattern(struct coder *coder,unsigned int pattern, int nbits, unsigned char **output); - -void Ptngc_writebits(struct coder *coder,unsigned int value,int nbits, unsigned char **output_ptr); -void Ptngc_write32bits(struct coder *coder,unsigned int value,int nbits, unsigned char **output_ptr); -void Ptngc_writemanybits(struct coder *coder,unsigned char *value,int nbits, unsigned char **output_ptr); +struct coder DECLSPECDLLEXPORT *Ptngc_coder_init(void); +void DECLSPECDLLEXPORT Ptngc_coder_deinit(struct coder *coder); +unsigned char DECLSPECDLLEXPORT *Ptngc_pack_array(struct coder *coder,int *input, int *length, int coding, int coding_parameter, int natoms, int speed); +int DECLSPECDLLEXPORT Ptngc_unpack_array(struct coder *coder,unsigned char *packed,int *output, int length, int coding, int coding_parameter, int natoms); +unsigned char DECLSPECDLLEXPORT *Ptngc_pack_array_xtc2(struct coder *coder,int *input, int *length); +int DECLSPECDLLEXPORT Ptngc_unpack_array_xtc2(struct coder *coder,unsigned char *packed,int *output, int length); +unsigned char DECLSPECDLLEXPORT *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int speed); +int DECLSPECDLLEXPORT Ptngc_unpack_array_xtc3(unsigned char *packed,int *output, int length, int natoms); + +void DECLSPECDLLEXPORT Ptngc_out8bits(struct coder *coder, unsigned char **output); +void DECLSPECDLLEXPORT Ptngc_pack_flush(struct coder *coder,unsigned char **output); +void DECLSPECDLLEXPORT Ptngc_write_pattern(struct coder *coder,unsigned int pattern, int nbits, unsigned char **output); + +void DECLSPECDLLEXPORT Ptngc_writebits(struct coder *coder,unsigned int value,int nbits, unsigned char **output_ptr); +void DECLSPECDLLEXPORT Ptngc_write32bits(struct coder *coder,unsigned int value,int nbits, unsigned char **output_ptr); +void DECLSPECDLLEXPORT Ptngc_writemanybits(struct coder *coder,unsigned char *value,int nbits, unsigned char **output_ptr); #endif diff --git a/src/compression/bwlzh.c b/src/compression/bwlzh.c index 783a253..c9ab553 100644 --- a/src/compression/bwlzh.c +++ b/src/compression/bwlzh.c @@ -66,15 +66,15 @@ static void printvals(char *name, unsigned int *vals, int nvals) fprintf(stderr,"\n"); } #endif - + } #endif static void bwlzh_compress_gen(unsigned int *vals, int nvals, - unsigned char *output, int *output_len, - int enable_lz77, - int verbose) + unsigned char *output, int *output_len, + int enable_lz77, + int verbose) { unsigned int *vals16; int nvals16; @@ -119,7 +119,7 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, lens=tmpmem+max_vals_per_block*15; #ifdef PARTIAL_MTF3 mtf3=warnmalloc(max_vals_per_block*3*3*sizeof *mtf3); /* 3 due to expansion of 32 bit to 16 bit, 3 due to up to 3 bytes - per 16 value. */ + per 16 value. */ #endif if (verbose) { @@ -131,20 +131,20 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, output[outdata++]=(((unsigned int)nvals)>>8)&0xFFU; output[outdata++]=(((unsigned int)nvals)>>16)&0xFFU; output[outdata++]=(((unsigned int)nvals)>>24)&0xFFU; - + valsleft=nvals; valstart=0; while (valsleft) { int reducealgo=1; /* Reduce algo is LZ77. */ if (!enable_lz77) - reducealgo=0; + reducealgo=0; thisvals=valsleft; if (thisvals>max_vals_per_block) - thisvals=max_vals_per_block; + thisvals=max_vals_per_block; valsleft-=thisvals; if (verbose) - fprintf(stderr,"Creating vals16 block from %d values.\n",thisvals); + fprintf(stderr,"Creating vals16 block from %d values.\n",thisvals); #ifdef SHOWIT printvals("vals",vals+valstart,thisvals); @@ -160,15 +160,15 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, #ifdef SHOWIT printvals("vals16",vals16,nvals16); #endif - + if (verbose) - { - fprintf(stderr,"Resulting vals16 values: %d\n",nvals16); - } + { + fprintf(stderr,"Resulting vals16 values: %d\n",nvals16); + } if (verbose) - { - fprintf(stderr,"BWT\n"); - } + { + fprintf(stderr,"BWT\n"); + } Ptngc_comp_to_bwt(vals16,nvals16,bwt,&bwt_index); #ifdef SHOWIT @@ -181,7 +181,7 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, output[outdata++]=(((unsigned int)thisvals)>>8)&0xFFU; output[outdata++]=(((unsigned int)thisvals)>>16)&0xFFU; output[outdata++]=(((unsigned int)thisvals)>>24)&0xFFU; - + /* Store the number of nvals16 values in this block. */ output[outdata++]=((unsigned int)nvals16)&0xFFU; output[outdata++]=(((unsigned int)nvals16)>>8)&0xFFU; @@ -193,261 +193,261 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, output[outdata++]=(((unsigned int)bwt_index)>>8)&0xFFU; output[outdata++]=(((unsigned int)bwt_index)>>16)&0xFFU; output[outdata++]=(((unsigned int)bwt_index)>>24)&0xFFU; - + if (verbose) - fprintf(stderr,"MTF\n"); + fprintf(stderr,"MTF\n"); #ifdef PARTIAL_MTF3 Ptngc_comp_conv_to_mtf_partial3(bwt,nvals16, - mtf3); + mtf3); for (imtfinner=0; imtfinner<3; imtfinner++) - { - int i; - if (verbose) - fprintf(stderr,"Doing partial MTF: %d\n",imtfinner); - for (i=0; i>=1; - } - coarse[numbits-1]+=thist[jj]; - } + if (noffsets) + { + unsigned int thist[0x20004]; + unsigned int coarse[17]={0,}; + int jj; + Ptngc_comp_make_dict_hist(lens,nlens,dict,&ndict,thist); + for (jj=0; jj>=1; + } + coarse[numbits-1]+=thist[jj]; + } #if 1 - for (jj=0; jj>8)&0xFFU; - output[outdata++]=(((unsigned int)nrle)>>16)&0xFFU; - output[outdata++]=(((unsigned int)nrle)>>24)&0xFFU; - - /* Store the size of the huffman block. */ - output[outdata++]=((unsigned int)bwlzhhufflen)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>8)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; - - /* Store the huffman block. */ - memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); - outdata+=bwlzhhufflen; - - if (reducealgo==1) - { - /* Store the number of values in this block. */ - output[outdata++]=((unsigned int)noffsets)&0xFFU; - output[outdata++]=(((unsigned int)noffsets)>>8)&0xFFU; - output[outdata++]=(((unsigned int)noffsets)>>16)&0xFFU; - output[outdata++]=(((unsigned int)noffsets)>>24)&0xFFU; - - if (noffsets>0) - { - if (verbose) - fprintf(stderr,"Huffman for offsets\n"); - - huffalgo=-1; - Ptngc_comp_huff_compress_verbose(offsets,noffsets,bwlzhhuff,&bwlzhhufflen,&huffdatalen,nhufflen,&huffalgo,1); - if (verbose) - { - int i; - fprintf(stderr,"Huffman data length is %d B.\n",huffdatalen); - for (i=0; i>8)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; - - /* Store the huffman block. */ - memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); - outdata+=bwlzhhufflen; - } - else - { - int i; - output[outdata++]=1; - for (i=0; i>8)&0xFFU; - } - if (verbose) - fprintf(stderr,"Store raw offsets: %d B\n",noffsets*2); - } - } + { + int i; + fprintf(stderr,"Huffman\n"); + for (i=0; i>8)&0xFFU; + output[outdata++]=(((unsigned int)nrle)>>16)&0xFFU; + output[outdata++]=(((unsigned int)nrle)>>24)&0xFFU; + + /* Store the size of the huffman block. */ + output[outdata++]=((unsigned int)bwlzhhufflen)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>8)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; + + /* Store the huffman block. */ + memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); + outdata+=bwlzhhufflen; + + if (reducealgo==1) + { + /* Store the number of values in this block. */ + output[outdata++]=((unsigned int)noffsets)&0xFFU; + output[outdata++]=(((unsigned int)noffsets)>>8)&0xFFU; + output[outdata++]=(((unsigned int)noffsets)>>16)&0xFFU; + output[outdata++]=(((unsigned int)noffsets)>>24)&0xFFU; + + if (noffsets>0) + { + if (verbose) + fprintf(stderr,"Huffman for offsets\n"); + + huffalgo=-1; + Ptngc_comp_huff_compress_verbose(offsets,noffsets,bwlzhhuff,&bwlzhhufflen,&huffdatalen,nhufflen,&huffalgo,1); + if (verbose) + { + int i; + fprintf(stderr,"Huffman data length is %d B.\n",huffdatalen); + for (i=0; i>8)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; + + /* Store the huffman block. */ + memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); + outdata+=bwlzhhufflen; + } + else + { + int i; + output[outdata++]=1; + for (i=0; i>8)&0xFFU; + } + if (verbose) + fprintf(stderr,"Store raw offsets: %d B\n",noffsets*2); + } + } #if 0 - { - int i,ndict; - FILE *f=fopen("len.dict","w"); - Ptngc_comp_make_dict_hist(lens,nlens,dict,&ndict,hist); - for (i=0; i>8)&0xFFU; - output[outdata++]=(((unsigned int)nlens)>>16)&0xFFU; - output[outdata++]=(((unsigned int)nlens)>>24)&0xFFU; - - /* Store the size of the huffman block. */ - output[outdata++]=((unsigned int)bwlzhhufflen)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>8)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; - output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; - - /* Store the huffman block. */ - memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); - outdata+=bwlzhhufflen; - } + { + int i,ndict; + FILE *f=fopen("len.dict","w"); + Ptngc_comp_make_dict_hist(lens,nlens,dict,&ndict,hist); + for (i=0; i>8)&0xFFU; + output[outdata++]=(((unsigned int)nlens)>>16)&0xFFU; + output[outdata++]=(((unsigned int)nlens)>>24)&0xFFU; + + /* Store the size of the huffman block. */ + output[outdata++]=((unsigned int)bwlzhhufflen)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>8)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>16)&0xFFU; + output[outdata++]=(((unsigned int)bwlzhhufflen)>>24)&0xFFU; + + /* Store the huffman block. */ + memcpy(output+outdata,bwlzhhuff,bwlzhhufflen); + outdata+=bwlzhhufflen; + } #ifdef PARTIAL_MTF3 - } + } #endif } @@ -463,39 +463,39 @@ static void bwlzh_compress_gen(unsigned int *vals, int nvals, void DECLSPECDLLEXPORT bwlzh_compress(unsigned int *vals, int nvals, - unsigned char *output, int *output_len) + unsigned char *output, int *output_len) { bwlzh_compress_gen(vals,nvals,output,output_len,1,0); } void DECLSPECDLLEXPORT bwlzh_compress_verbose(unsigned int *vals, int nvals, - unsigned char *output, int *output_len) + unsigned char *output, int *output_len) { bwlzh_compress_gen(vals,nvals,output,output_len,1,1); } void DECLSPECDLLEXPORT bwlzh_compress_no_lz77(unsigned int *vals, int nvals, - unsigned char *output, int *output_len) + unsigned char *output, int *output_len) { bwlzh_compress_gen(vals,nvals,output,output_len,0,0); } void DECLSPECDLLEXPORT bwlzh_compress_no_lz77_verbose(unsigned int *vals, int nvals, - unsigned char *output, int *output_len) + unsigned char *output, int *output_len) { bwlzh_compress_gen(vals,nvals,output,output_len,0,1); } static void bwlzh_decompress_gen(unsigned char *input, int nvals, - unsigned int *vals, - int verbose) + unsigned int *vals, + int verbose) { unsigned int *vals16; int nvals16; int bwt_index; - unsigned int *bwt=NULL; + unsigned int *bwt=NULL; unsigned int *mtf=NULL; #ifdef PARTIAL_MTF3 unsigned char *mtf3=NULL; @@ -532,7 +532,7 @@ static void bwlzh_decompress_gen(unsigned char *input, int nvals, lens=tmpmem+max_vals_per_block*15; #ifdef PARTIAL_MTF3 mtf3=warnmalloc(max_vals_per_block*3*3*sizeof *mtf3); /* 3 due to expansion of 32 bit to 16 bit, 3 due to up to 3 bytes - per 16 value. */ + per 16 value. */ #endif if (verbose) @@ -542,9 +542,9 @@ static void bwlzh_decompress_gen(unsigned char *input, int nvals, /* Read the number of real values in the whole block. */ nvalsfile=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); + (((unsigned int)input[inpdata+1])<<8) | + (((unsigned int)input[inpdata+2])<<16) | + (((unsigned int)input[inpdata+3])<<24)); inpdata+=4; if (nvalsfile!=nvals) @@ -552,7 +552,7 @@ static void bwlzh_decompress_gen(unsigned char *input, int nvals, fprintf(stderr,"BWLZH: The number of values found in the file is different from the number of values expected.\n"); exit(EXIT_FAILURE); } - + valsleft=nvals; valstart=0; while (valsleft) @@ -561,182 +561,182 @@ static void bwlzh_decompress_gen(unsigned char *input, int nvals, int reducealgo; /* Read the number of real values in this block. */ thisvals=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); + (((unsigned int)input[inpdata+1])<<8) | + (((unsigned int)input[inpdata+2])<<16) | + (((unsigned int)input[inpdata+3])<<24)); inpdata+=4; - + valsleft-=thisvals; /* Read the number of nvals16 values in this block. */ nvals16=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); + (((unsigned int)input[inpdata+1])<<8) | + (((unsigned int)input[inpdata+2])<<16) | + (((unsigned int)input[inpdata+3])<<24)); inpdata+=4; /* Read the BWT index. */ bwt_index=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); + (((unsigned int)input[inpdata+1])<<8) | + (((unsigned int)input[inpdata+2])<<16) | + (((unsigned int)input[inpdata+3])<<24)); inpdata+=4; if (thisvals>max_vals_per_block) - { - /* More memory must be allocated for decompression. */ - max_vals_per_block=thisvals; - if (verbose) - fprintf(stderr,"Allocating more memory: %d B\n",(int)(max_vals_per_block*15*sizeof *tmpmem)); - tmpmem=warnrealloc(tmpmem,max_vals_per_block*18*sizeof *tmpmem); - vals16=tmpmem; - bwt=tmpmem+max_vals_per_block*3; - mtf=tmpmem+max_vals_per_block*6; - rle=tmpmem+max_vals_per_block*9; - offsets=tmpmem+max_vals_per_block*12; - lens=tmpmem+max_vals_per_block*15; + { + /* More memory must be allocated for decompression. */ + max_vals_per_block=thisvals; + if (verbose) + fprintf(stderr,"Allocating more memory: %d B\n",(int)(max_vals_per_block*15*sizeof *tmpmem)); + tmpmem=warnrealloc(tmpmem,max_vals_per_block*18*sizeof *tmpmem); + vals16=tmpmem; + bwt=tmpmem+max_vals_per_block*3; + mtf=tmpmem+max_vals_per_block*6; + rle=tmpmem+max_vals_per_block*9; + offsets=tmpmem+max_vals_per_block*12; + lens=tmpmem+max_vals_per_block*15; #ifdef PARTIAL_MTF3 - mtf3=warnrealloc(mtf3,max_vals_per_block*3*3*sizeof *mtf3); /* 3 due to expansion of 32 bit to 16 bit, 3 due to up to 3 bytes - per 16 value. */ -#endif - } + mtf3=warnrealloc(mtf3,max_vals_per_block*3*3*sizeof *mtf3); /* 3 due to expansion of 32 bit to 16 bit, 3 due to up to 3 bytes + per 16 value. */ +#endif + } #ifdef PARTIAL_MTF3 for (imtfinner=0; imtfinner<3; imtfinner++) - { - int i; - if (verbose) - fprintf(stderr,"Doing partial MTF: %d\n",imtfinner); -#endif - - reducealgo=(int)input[inpdata]; - inpdata++; - - /* Read the number of huffman values in this block. */ - nrle=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); - inpdata+=4; - - /* Read the size of the huffman block. */ - bwlzhhufflen=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); - inpdata+=4; - - if (verbose) - fprintf(stderr,"Decompressing huffman block of length %d.\n",bwlzhhufflen); - /* Decompress the huffman block. */ - Ptngc_comp_huff_decompress(input+inpdata,bwlzhhufflen,rle); - inpdata+=bwlzhhufflen; - - if (reducealgo==1) /* LZ77 */ - { - int offstore; - /* Read the number of huffman values in this block. */ - noffsets=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); - inpdata+=4; - - if (noffsets>0) - { - /* How are the offsets stored? */ - offstore=(int)input[inpdata++]; - if (offstore==0) - { - /* Read the size of the huffman block. */ - bwlzhhufflen=(int)(((unsigned int)input[inpdata]) | - (((unsigned int)input[inpdata+1])<<8) | - (((unsigned int)input[inpdata+2])<<16) | - (((unsigned int)input[inpdata+3])<<24)); - inpdata+=4; - - if (verbose) - fprintf(stderr,"Decompressing offset huffman block.\n"); - - /* Decompress the huffman block. */ - Ptngc_comp_huff_decompress(input+inpdata,bwlzhhufflen,offsets); - inpdata+=bwlzhhufflen; - } - else - { - int i; - if (verbose) - fprintf(stderr,"Reading offset block.\n"); - for (i=0; i0) + { + /* How are the offsets stored? */ + offstore=(int)input[inpdata++]; + if (offstore==0) + { + /* Read the size of the huffman block. */ + bwlzhhufflen=(int)(((unsigned int)input[inpdata]) | + (((unsigned int)input[inpdata+1])<<8) | + (((unsigned int)input[inpdata+2])<<16) | + (((unsigned int)input[inpdata+3])<<24)); + inpdata+=4; + + if (verbose) + fprintf(stderr,"Decompressing offset huffman block.\n"); + + /* Decompress the huffman block. */ + Ptngc_comp_huff_decompress(input+inpdata,bwlzhhufflen,offsets); + inpdata+=bwlzhhufflen; + } + else + { + int i; + if (verbose) + fprintf(stderr,"Reading offset block.\n"); + for (i=0; i>8); int k1=(int)(nrepeat[i1]&0xFFU); int repeat2=(int)(nrepeat[i2]>>8); int k2=(int)(nrepeat[i2]&0xFFU); if ((repeat1>1) && (repeat2>1) && (k1==k2)) - { - int maxskip=0; - /* Yes. Compare the repeating patterns. */ - for (j=0; jv2) - return 1; - } - /* The repeating patterns are equal. Skip as far as we can - before continuing. */ - maxskip=repeat1; - if (repeat2v2) + return 1; + } + /* The repeating patterns are equal. Skip as far as we can + before continuing. */ + maxskip=repeat1; + if (repeat2vals[i2]) - return 1; - i1++; - if (i1>=nvals) - i1=0; - i2++; - if (i2>=nvals) - i2=0; - } + { + if (vals[i1]vals[i2]) + return 1; + i1++; + if (i1>=nvals) + i1=0; + i2++; + if (i2>=nvals) + i2=0; + } } return 0; } void Ptngc_bwt_merge_sort_inner(int *indices, int nvals,unsigned int *vals, - int start, int end, - unsigned int *nrepeat, - int *workarray) + int start, int end, + unsigned int *nrepeat, + int *workarray) { int middle; if ((end-start)>1) @@ -90,60 +90,60 @@ void Ptngc_bwt_merge_sort_inner(int *indices, int nvals,unsigned int *vals, printf("For start %d end %d obtained new middle: %d\n",start,end,middle); #endif Ptngc_bwt_merge_sort_inner(indices,nvals,vals, - start,middle, - nrepeat, - workarray); + start,middle, + nrepeat, + workarray); Ptngc_bwt_merge_sort_inner(indices,nvals,vals, - middle,end, - nrepeat, - workarray); + middle,end, + nrepeat, + workarray); #if 0 printf("For start %d end %d Before merge: Comparing element %d with %d\n",start,end,middle-1,middle); #endif if (compare_index(indices[middle-1],indices[middle],nvals,vals,nrepeat)>0) - { - /* Merge to work array. */ - int i, n=end-start; - int ileft=start; - int iright=middle; - for (i=0; i0) - { - workarray[i]=indices[iright]; - iright++; - } - else - { - workarray[i]=indices[ileft]; - ileft++; - } - } - } - /* Copy result back. */ - memcpy(indices+start,workarray,(end-start)*sizeof(int)); - } + if (compare_index(indices[ileft],indices[iright],nvals,vals,nrepeat)>0) + { + workarray[i]=indices[iright]; + iright++; + } + else + { + workarray[i]=indices[ileft]; + ileft++; + } + } + } + /* Copy result back. */ + memcpy(indices+start,workarray,(end-start)*sizeof(int)); + } } } /* Burrows-Wheeler transform. */ void Ptngc_comp_to_bwt(unsigned int *vals, int nvals, - unsigned int *output, int *index) + unsigned int *output, int *index) { int i; int *indices=warnmalloc(2*nvals*sizeof *indices); @@ -172,94 +172,94 @@ void Ptngc_comp_to_bwt(unsigned int *vals, int nvals, for (i=0; i=1; k--) - { - try_next_k: - if (k>=1) - { + { + int maxrepeat=nvals*2; + int j,k,m; + int good_j=-1, good_k=0; + int kmax=16; + /* Track repeating patterns. + k=1 corresponds to AAAAA... + k=2 corresponds to ABABAB... + k=3 corresponds to ABCABCABCABC... + k=4 corresponds to ABCDABCDABCD... + etc. */ + for (k=kmax; k>=1; k--) + { + try_next_k: + if (k>=1) + { #ifdef SHOWIT - printf("Trying k=%d at i=%d\n",k,i); + printf("Trying k=%d at i=%d\n",k,i); #endif - for (j=k; jmaxrepeat) - new_j=j; - if ((new_j>good_j) || ((new_j==good_j) && (kmaxrepeat) + new_j=j; + if ((new_j>good_j) || ((new_j==good_j) && (knvals) - repeat=nvals; - nrepeat[i+m]=((unsigned int) (good_k)) | (((unsigned int) (repeat))<<8); - } - /* If no repetition was found for this value signal that here. */ - if (!nrepeat[i]) - nrepeat[i+m]=257U; /* This is 1<<8 | 1 */ - } + goto try_next_k; + } + } + } + } + /* From good_j and good_k we know the repeat for a large + number of strings. The very last repeat length should not + be assigned, since it can be much longer if a new test is + done. */ + for (m=0; (m+good_knvals) + repeat=nvals; + nrepeat[i+m]=((unsigned int) (good_k)) | (((unsigned int) (repeat))<<8); + } + /* If no repetition was found for this value signal that here. */ + if (!nrepeat[i]) + nrepeat[i+m]=257U; /* This is 1<<8 | 1 */ + } } #ifdef SHOWIT for (i=0; ipack_temporary_bits=0; return coder_inst; } -void Ptngc_coder_deinit(struct coder *coder_inst) +void DECLSPECDLLEXPORT Ptngc_coder_deinit(struct coder *coder_inst) { free(coder_inst); } -TNG_INLINE void Ptngc_out8bits(struct coder *coder_inst, unsigned char **output) +TNG_INLINE void DECLSPECDLLEXPORT Ptngc_out8bits(struct coder *coder_inst, unsigned char **output) { int pack_temporary_bits=coder_inst->pack_temporary_bits; unsigned int pack_temporary=coder_inst->pack_temporary; @@ -61,7 +61,7 @@ TNG_INLINE void Ptngc_out8bits(struct coder *coder_inst, unsigned char **output) coder_inst->pack_temporary=pack_temporary; } -void Ptngc_write_pattern(struct coder *coder_inst, unsigned int pattern, +void DECLSPECDLLEXPORT Ptngc_write_pattern(struct coder *coder_inst, unsigned int pattern, int nbits, unsigned char **output) { unsigned int mask1,mask2; @@ -71,17 +71,17 @@ void Ptngc_write_pattern(struct coder *coder_inst, unsigned int pattern, coder_inst->pack_temporary_bits+=nbits; while (nbits) { - if (pattern & mask1) - coder_inst->pack_temporary|=mask2; - nbits--; - mask1<<=1; - mask2>>=1; + if (pattern & mask1) + coder_inst->pack_temporary|=mask2; + nbits--; + mask1<<=1; + mask2>>=1; } Ptngc_out8bits(coder_inst,output); } /* Write up to 24 bits */ -TNG_INLINE void Ptngc_writebits(struct coder *coder_inst, +TNG_INLINE void DECLSPECDLLEXPORT Ptngc_writebits(struct coder *coder_inst, unsigned int value, int nbits, unsigned char **output_ptr) { @@ -93,7 +93,7 @@ TNG_INLINE void Ptngc_writebits(struct coder *coder_inst, } /* Write up to 32 bits */ -void Ptngc_write32bits(struct coder *coder_inst,unsigned int value, +void DECLSPECDLLEXPORT Ptngc_write32bits(struct coder *coder_inst,unsigned int value, int nbits, unsigned char **output_ptr) { unsigned int mask; @@ -116,15 +116,15 @@ void Ptngc_write32bits(struct coder *coder_inst,unsigned int value, } /* Write "arbitrary" number of bits */ -void Ptngc_writemanybits(struct coder *coder_inst, unsigned char *value, +void DECLSPECDLLEXPORT Ptngc_writemanybits(struct coder *coder_inst, unsigned char *value, int nbits, unsigned char **output_ptr) { int vptr=0; while (nbits>=24) { unsigned int v=((((unsigned int)value[vptr])<<16)| - (((unsigned int)value[vptr+1])<<8)| - (((unsigned int)value[vptr+2]))); + (((unsigned int)value[vptr+1])<<8)| + (((unsigned int)value[vptr+2]))); Ptngc_writebits(coder_inst,v,24,output_ptr); vptr+=3; nbits-=24; @@ -151,8 +151,8 @@ static int write_stop_bit_code(struct coder *coder_inst, unsigned int s, s>>=coding_parameter; if (s) { - this|=1U; - coder_inst->stat_overflow++; + this|=1U; + coder_inst->stat_overflow++; } coder_inst->pack_temporary<<=(coding_parameter+1); coder_inst->pack_temporary_bits+=coding_parameter+1; @@ -160,9 +160,9 @@ static int write_stop_bit_code(struct coder *coder_inst, unsigned int s, Ptngc_out8bits(coder_inst,output); if (s) { - coding_parameter>>=1; - if (coding_parameter<1) - coding_parameter=1; + coding_parameter>>=1; + if (coding_parameter<1) + coding_parameter=1; } } while (s); coder_inst->stat_numval++; @@ -175,15 +175,15 @@ static int pack_stopbits_item(struct coder *coder_inst,int item, /* Find this symbol in table. */ int s=0; if (item>0) - s=1+(item-1)*2; + s=1+(item-1)*2; else if (item<0) - s=2+(-item-1)*2; + s=2+(-item-1)*2; return write_stop_bit_code(coder_inst,s,coding_parameter,output); } static int pack_triplet(struct coder *coder_inst, unsigned int *s, unsigned char **output, int coding_parameter, - unsigned int max_base, int maxbits) + unsigned int max_base, int maxbits) { /* Determine base for this triplet. */ unsigned int min_base=1U<=this_base) { - this_base*=2; - jbase++; + this_base*=2; + jbase++; } bits_per_value=coding_parameter+jbase; if (jbase>=3) { if (this_base>max_base) - return 1; + return 1; bits_per_value=maxbits; jbase=3; } @@ -215,14 +215,14 @@ static int pack_triplet(struct coder *coder_inst, unsigned int *s, return 0; } -void Ptngc_pack_flush(struct coder *coder_inst,unsigned char **output) +void DECLSPECDLLEXPORT Ptngc_pack_flush(struct coder *coder_inst,unsigned char **output) { /* Zero-fill just enough. */ if (coder_inst->pack_temporary_bits>0) Ptngc_write_pattern(coder_inst,0,8-coder_inst->pack_temporary_bits,output); } -unsigned char *Ptngc_pack_array(struct coder *coder_inst, +unsigned char DECLSPECDLLEXPORT *Ptngc_pack_array(struct coder *coder_inst, int *input, int *length, int coding, int coding_parameter, int natoms, int speed) { @@ -235,25 +235,25 @@ unsigned char *Ptngc_pack_array(struct coder *coder_inst, int cnt=0; int most_negative=2147483647; for (i=0; i>8)&0xFFU; output[2]=(((unsigned int)most_negative)>>16)&0xFFU; output[3]=(((unsigned int)most_negative)>>24)&0xFFU; for (i=0; i=5) - bwlzh_compress(pval,n,output+4,length); + bwlzh_compress(pval,n,output+4,length); else - bwlzh_compress_no_lz77(pval,n,output+4,length); + bwlzh_compress_no_lz77(pval,n,output+4,length); (*length)+=4; free(pval); return output; @@ -275,64 +275,64 @@ unsigned char *Ptngc_pack_array(struct coder *coder_inst, output=warnmalloc(8* *length*sizeof *output); output_ptr=output; if ((coding==TNG_COMPRESS_ALGO_TRIPLET) || - (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA) || - (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE)) - { - /* Pack triplets. */ - int ntriplets=*length/3; - /* Determine max base and maxbits */ - unsigned int max_base=1U<0) - s=1+(item-1)*2; - else if (item<0) - s=2+(-item-1)*2; - if (s>intmax) - intmax=s; - } - /* Store intmax */ - coder_inst->pack_temporary_bits=32; - coder_inst->pack_temporary=intmax; - Ptngc_out8bits(coder_inst,&output_ptr); - while (intmax>=max_base) - { - max_base*=2; - maxbits++; - } - for (i=0; i0) - s[j]=1+(item-1)*2; - else if (item<0) - s[j]=2+(-item-1)*2; - } - if (pack_triplet(coder_inst, s, &output_ptr, + (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA) || + (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE)) + { + /* Pack triplets. */ + int ntriplets=*length/3; + /* Determine max base and maxbits */ + unsigned int max_base=1U<0) + s=1+(item-1)*2; + else if (item<0) + s=2+(-item-1)*2; + if (s>intmax) + intmax=s; + } + /* Store intmax */ + coder_inst->pack_temporary_bits=32; + coder_inst->pack_temporary=intmax; + Ptngc_out8bits(coder_inst,&output_ptr); + while (intmax>=max_base) + { + max_base*=2; + maxbits++; + } + for (i=0; i0) + s[j]=1+(item-1)*2; + else if (item<0) + s[j]=2+(-item-1)*2; + } + if (pack_triplet(coder_inst, s, &output_ptr, coding_parameter, max_base,maxbits)) - { - free(output); - return NULL; - } - } - } + { + free(output); + return NULL; + } + } + } else - for (i=0; i<*length; i++) - if (pack_stopbits_item(coder_inst,input[i],&output_ptr,coding_parameter)) - { - free(output); - return NULL; - } + for (i=0; i<*length; i++) + if (pack_stopbits_item(coder_inst,input[i],&output_ptr,coding_parameter)) + { + free(output); + return NULL; + } Ptngc_pack_flush(coder_inst,&output_ptr); output_length=(int)(output_ptr-output); *length=output_length; @@ -357,39 +357,39 @@ static int unpack_array_stop_bits(struct coder *coder_inst, unsigned int insert_mask=1U<<(numbits-1); int inserted_bits=numbits; do { - for (j=0; j>=1; - extract_mask>>=1; - if (!extract_mask) - { - extract_mask=0x80; - ptr++; - } - } - /* Check stop bit */ - bit=*ptr & extract_mask; - extract_mask>>=1; - if (!extract_mask) - { - extract_mask=0x80; - ptr++; - } - if (bit) - { - numbits>>=1; - if (numbits<1) - numbits=1; - inserted_bits+=numbits; - insert_mask=1U<<(inserted_bits-1); - } + for (j=0; j>=1; + extract_mask>>=1; + if (!extract_mask) + { + extract_mask=0x80; + ptr++; + } + } + /* Check stop bit */ + bit=*ptr & extract_mask; + extract_mask>>=1; + if (!extract_mask) + { + extract_mask=0x80; + ptr++; + } + if (bit) + { + numbits>>=1; + if (numbits<1) + numbits=1; + inserted_bits+=numbits; + insert_mask=1U<<(inserted_bits-1); + } } while (bit); s=(pattern+1)/2; if ((pattern%2)==0) - s=-s; + s=-s; output[i]=s; } return 0; @@ -426,45 +426,45 @@ static int unpack_array_triplet(struct coder *coder_inst, unsigned int numbits; unsigned int bit; for (j=0; j<2; j++) - { - bit=*ptr & extract_mask; - jbase<<=1; - if (bit) - jbase|=1U; - extract_mask>>=1; - if (!extract_mask) - { - extract_mask=0x80; - ptr++; - } - } + { + bit=*ptr & extract_mask; + jbase<<=1; + if (bit) + jbase|=1U; + extract_mask>>=1; + if (!extract_mask) + { + extract_mask=0x80; + ptr++; + } + } if (jbase==3) - numbits=maxbits; + numbits=maxbits; else - numbits=coding_parameter+jbase; + numbits=coding_parameter+jbase; for (j=0; j<3; j++) - { - int s; - unsigned int jbit; - unsigned int pattern=0; - for (jbit=0; jbit>=1; - if (!extract_mask) - { - extract_mask=0x80; - ptr++; - } - } - s=(pattern+1)/2; - if ((pattern%2)==0) - s=-s; - output[i*3+j]=s; - } + { + int s; + unsigned int jbit; + unsigned int pattern=0; + for (jbit=0; jbit>=1; + if (!extract_mask) + { + extract_mask=0x80; + ptr++; + } + } + s=(pattern+1)/2; + if ((pattern%2)==0) + s=-s; + output[i*3+j]=s; + } } return 0; } @@ -478,23 +478,23 @@ static int unpack_array_bwlzh(struct coder *coder_inst, int nframes=n/natoms/3; int cnt=0; int most_negative=(int)(((unsigned int)packed[0]) | - (((unsigned int)packed[1])<<8) | - (((unsigned int)packed[2])<<16) | - (((unsigned int)packed[3])<<24)); + (((unsigned int)packed[1])<<8) | + (((unsigned int)packed[2])<<16) | + (((unsigned int)packed[3])<<24)); (void) coder_inst; bwlzh_decompress(packed+4,length,pval); for (i=0; ileaf.idict].code=(code<<1)|htree->leaf.bit; #if 0 printf("I am a leaf: %d %d\n", - codelength[htree->leaf.idict].length, - codelength[htree->leaf.idict].code); + codelength[htree->leaf.idict].length, + codelength[htree->leaf.idict].code); #endif } else { if (!top) - { - code<<=1; - code|=htree->node.bit; - length++; - } + { + code<<=1; + code|=htree->node.bit; + length++; + } #if 0 printf("I am a node length: %d\n",length); printf("I am a node code: %d\n",code); @@ -109,14 +109,14 @@ static void free_nodes(union htree_nodeleaf *htree, int top) if (htree->nodeleaf==htree_leaf) { if (!top) - free(htree); + free(htree); } else { free_nodes(htree->node.n1,0); free_nodes(htree->node.n2,0); if (!top) - free(htree); + free(htree); } } @@ -174,12 +174,12 @@ static unsigned int readbits(int length, unsigned char **input, int *bitptr) *bitptr=(*bitptr)+1; extract_mask>>=1; if (!extract_mask) - { - extract_mask=0x80U; - *input=(*input)+1; - *bitptr=0; - thisval=**input; - } + { + extract_mask=0x80U; + *input=(*input)+1; + *bitptr=0; + thisval=**input; + } } return val; } @@ -219,14 +219,14 @@ static int comp_codes_value(const void *codeptr1, const void *codeptr2, const vo huffman_dict_unpacked array should be 131077 long (note five longer than 0x20000) */ void Ptngc_comp_conv_to_huffman(unsigned int *vals, int nvals, - unsigned int *dict, int ndict, - unsigned int *prob, - unsigned char *huffman, - int *huffman_len, - unsigned char *huffman_dict, - int *huffman_dictlen, - unsigned int *huffman_dict_unpacked, - int *huffman_dict_unpackedlen) + unsigned int *dict, int ndict, + unsigned int *prob, + unsigned char *huffman, + int *huffman_len, + unsigned char *huffman_dict, + int *huffman_dictlen, + unsigned int *huffman_dict_unpacked, + int *huffman_dict_unpackedlen) { int i; int nleft; @@ -239,121 +239,121 @@ void Ptngc_comp_conv_to_huffman(unsigned int *vals, int nvals, while (longcodes) { /* Create array of leafs (will be array of nodes/trees during - buildup of tree. */ + buildup of tree. */ htree=warnmalloc(ndict*sizeof *htree); codelength=warnmalloc(ndict*sizeof *codelength); bitptr=0; huffman_ptr=huffman; for (i=0; i1) - { - union htree_nodeleaf *n1=warnmalloc(sizeof *n1); - union htree_nodeleaf *n2=warnmalloc(sizeof *n2); - int new_place; - int p1,p2, new_prob; - *n1=htree[nleft-1]; - *n2=htree[nleft-2]; - if (n1->nodeleaf==htree_leaf) - { - p1=n1->leaf.prob; - n1->leaf.bit=0; - } - else - { - p1=n1->node.prob; - n1->node.bit=0; - } - if (n2->nodeleaf==htree_leaf) - { - p2=n2->leaf.prob; - n2->leaf.bit=1; - } - else - { - p2=n2->node.prob; - n2->node.bit=1; - } - nleft--; - /* Create a new node */ - htree[nleft-1].nodeleaf=htree_node; - htree[nleft-1].node.n1=n1; - htree[nleft-1].node.n2=n2; - new_prob=p1+p2; - htree[nleft-1].node.prob=new_prob; - /* Use insertion sort to place this in the correct place in the - array. */ - /* Where should it be inserted? */ - new_place=nleft; - while (new_place>0) - { - int pc; - if (htree[new_place-1].nodeleaf==htree_node) - pc=htree[new_place-1].node.prob; - else - pc=htree[new_place-1].leaf.prob; - if (new_prob1) + { + union htree_nodeleaf *n1=warnmalloc(sizeof *n1); + union htree_nodeleaf *n2=warnmalloc(sizeof *n2); + int new_place; + int p1,p2, new_prob; + *n1=htree[nleft-1]; + *n2=htree[nleft-2]; + if (n1->nodeleaf==htree_leaf) + { + p1=n1->leaf.prob; + n1->leaf.bit=0; + } + else + { + p1=n1->node.prob; + n1->node.bit=0; + } + if (n2->nodeleaf==htree_leaf) + { + p2=n2->leaf.prob; + n2->leaf.bit=1; + } + else + { + p2=n2->node.prob; + n2->node.bit=1; + } + nleft--; + /* Create a new node */ + htree[nleft-1].nodeleaf=htree_node; + htree[nleft-1].node.n1=n1; + htree[nleft-1].node.n2=n2; + new_prob=p1+p2; + htree[nleft-1].node.prob=new_prob; + /* Use insertion sort to place this in the correct place in the + array. */ + /* Where should it be inserted? */ + new_place=nleft; + while (new_place>0) + { + int pc; + if (htree[new_place-1].nodeleaf==htree_node) + pc=htree[new_place-1].node.prob; + else + pc=htree[new_place-1].leaf.prob; + if (new_probMAX_HUFFMAN_LEN) - longcodes=1; + if (codelength[i].length>MAX_HUFFMAN_LEN) + longcodes=1; /* If the codes are too long alter the probabilities. */ if (longcodes) - { - for (i=0; i>=1; - if (prob[i]==0) - prob[i]=1; - } + { + for (i=0; i>=1; + if (prob[i]==0) + prob[i]=1; + } - /* Free codelength. We will compute a new one. */ - free(codelength); - } + /* Free codelength. We will compute a new one. */ + free(codelength); + } } #if 0 { for (i=0; i=0; j--) - { - int bit=c&mask; - if (bit) - printf("1"); - else - printf("0"); - mask>>=1; - } - printf("\n"); - } + printf("%d %d\t\t %d %d ",codelength[i].dict,codelength[i].prob,codelength[i].length,codelength[i].code); + { + unsigned int c=codelength[i].code; + int j; + unsigned int mask=1<<(codelength[i].length-1); + for (j=codelength[i].length-1; j>=0; j--) + { + int bit=c&mask; + if (bit) + printf("1"); + else + printf("0"); + mask>>=1; + } + printf("\n"); + } } } #endif @@ -409,8 +409,8 @@ void Ptngc_comp_conv_to_huffman(unsigned int *vals, int nvals, { int r; for (r=0; r=0; j--) - { - int bit=c&mask; - if (bit) - printf("1"); - else - printf("0"); - mask>>=1; - } - printf("\n"); - } + printf("%d -\t\t %d %d ",codelength[i].dict,codelength[i].length,codelength[i].code); + { + unsigned int c=codelength[i].code; + int j; + unsigned int mask=1<<(codelength[i].length-1); + for (j=codelength[i].length-1; j>=0; j--) + { + int bit=c&mask; + if (bit) + printf("1"); + else + printf("0"); + mask>>=1; + } + printf("\n"); + } } } #endif @@ -565,17 +565,17 @@ void Ptngc_comp_conv_from_huffman(unsigned char *huffman, symbol=readbits(len,&huffman_ptr,&bitptr); j=0; while (symbol!=codelength[j].code) - { - int newlen; - j++; - newlen=codelength[j].length; - if (newlen!=len) - { - symbol<<=(newlen-len); - symbol|=readbits(newlen-len,&huffman_ptr,&bitptr); - len=newlen; - } - } + { + int newlen; + j++; + newlen=codelength[j].length; + if (newlen!=len) + { + symbol<<=(newlen-len); + symbol|=readbits(newlen-len,&huffman_ptr,&bitptr); + len=newlen; + } + } vals[i]=codelength[j].dict; } /* Free info about codes and length. */ diff --git a/src/compression/huffmem.c b/src/compression/huffmem.c index 2232b7e..dec6837 100644 --- a/src/compression/huffmem.c +++ b/src/compression/huffmem.c @@ -28,10 +28,10 @@ int Ptngc_comp_huff_buflen(int nvals) /* the value pointed to by chosen_algo should be sent as -1 for autodetect. */ void Ptngc_comp_huff_compress_verbose(unsigned int *vals, int nvals, - unsigned char *huffman, int *huffman_len, - int *huffdatalen, - int *huffman_lengths,int *chosen_algo, - int isvals16) + unsigned char *huffman, int *huffman_len, + int *huffdatalen, + int *huffman_lengths,int *chosen_algo, + int isvals16) { unsigned int *dict=warnmalloc(0x20005*sizeof *dict); unsigned int *hist=warnmalloc(0x20005*sizeof *hist); @@ -68,9 +68,9 @@ void Ptngc_comp_huff_compress_verbose(unsigned int *vals, int nvals, /* First compress the data using huffman coding (place it ready for output at 14 (code for algorithm+length etc.). */ Ptngc_comp_conv_to_huffman(vals,nvals,dict,ndict,hist, - huffman+14,&nhuff, - huffdict,&nhuffdict, - huffdictunpack,&nhuffdictunpack); + huffman+14,&nhuff, + huffdict,&nhuffdict, + huffdictunpack,&nhuffdictunpack); *huffdatalen=nhuff; /* Algorithm 0 stores the huffman dictionary directly (+ a code for @@ -83,30 +83,30 @@ void Ptngc_comp_huff_compress_verbose(unsigned int *vals, int nvals, Ptngc_comp_make_dict_hist(huffdictunpack,nhuffdictunpack,dict,&ndict1,hist); /* Pack huffman dictionary */ Ptngc_comp_conv_to_huffman(huffdictunpack,nhuffdictunpack, - dict,ndict1,hist, - huffman1,&nhuff1, - huffdict1,&nhuffdict1, - huffdictunpack1,&nhuffdictunpack1); + dict,ndict1,hist, + huffman1,&nhuff1, + huffdict1,&nhuffdict1, + huffdictunpack1,&nhuffdictunpack1); huffman_lengths[1]=nhuff+nhuff1+nhuffdict1+1*2+3*4+3+3+3+3+3; /* ... and rle + huffman coding ... (algorithm 2) Pack any repetetitive patterns. */ Ptngc_comp_conv_to_rle(huffdictunpack,nhuffdictunpack, - huffdictrle,&nhuffrle,1); + huffdictrle,&nhuffrle,1); /* Determine probabilities. */ Ptngc_comp_make_dict_hist(huffdictrle,nhuffrle,dict,&ndict2,hist); /* Pack huffman dictionary */ Ptngc_comp_conv_to_huffman(huffdictrle,nhuffrle, - dict,ndict2,hist, - huffman2,&nhuff2, - huffdict2,&nhuffdict2, - huffdictunpack2,&nhuffdictunpack2); + dict,ndict2,hist, + huffman2,&nhuff2, + huffdict2,&nhuffdict2, + huffdictunpack2,&nhuffdictunpack2); huffman_lengths[2]=nhuff+nhuff2+nhuffdict2+1*2+3*4+3+3+3+3+3+3; /* Choose the best algorithm and output the data. */ if ((*chosen_algo==0) || ((*chosen_algo==-1) && - (((huffman_lengths[0]>8)&0xFFU; huffman[19+nhuff]=(((unsigned int)ndict)>>16)&0xFFU; for (i=0; i>8)&0xFFU; huffman[28+nhuff]=(((unsigned int)ndict1)>>16)&0xFFU; for (i=0; i>8)&0xFFU; huffman[31+nhuff]=(((unsigned int)ndict2)>>16)&0xFFU; for (i=0; i=1; k--) - { - try_next_k: - if (k>=1) - { - for (j=k; jmaxrepeat) - new_j=j; - if ((new_j>good_j) || ((new_j==good_j) && (knvals) - repeat=nvals; - nrepeat[i+m]=((unsigned int) (good_k)) | (((unsigned int) (repeat))<<8); - } - /* If no repetition was found for this value signal that here. */ - if (!nrepeat[i]) - nrepeat[i+m]=257U; /* This is 1<<8 | 1 */ - } + { + int maxrepeat=nvals*2; + int j,k,m; + int good_j=-1, good_k=0; + int kmax=16; + /* Track repeating patterns. + k=1 corresponds to AAAAA... + k=2 corresponds to ABABAB... + k=3 corresponds to ABCABCABCABC... + k=4 corresponds to ABCDABCDABCD... + etc. */ + for (k=kmax; k>=1; k--) + { + try_next_k: + if (k>=1) + { + for (j=k; jmaxrepeat) + new_j=j; + if ((new_j>good_j) || ((new_j==good_j) && (knvals) + repeat=nvals; + nrepeat[i+m]=((unsigned int) (good_k)) | (((unsigned int) (repeat))<<8); + } + /* If no repetition was found for this value signal that here. */ + if (!nrepeat[i]) + nrepeat[i+m]=257U; /* This is 1<<8 | 1 */ + } } /* Sort cyclic shift matrix. */ @@ -155,19 +155,19 @@ static void add_circular(int *previous,int v, int i) { previous[(NUM_PREVIOUS+3)*v]++; if (previous[(NUM_PREVIOUS+3)*v]>NUM_PREVIOUS) - previous[(NUM_PREVIOUS+3)*v]=NUM_PREVIOUS; + previous[(NUM_PREVIOUS+3)*v]=NUM_PREVIOUS; previous[(NUM_PREVIOUS+3)*v+3+previous[(NUM_PREVIOUS+3)*v+1]]=i; previous[(NUM_PREVIOUS+3)*v+1]++; if (previous[(NUM_PREVIOUS+3)*v+1]>=NUM_PREVIOUS) - previous[(NUM_PREVIOUS+3)*v+1]=0; + previous[(NUM_PREVIOUS+3)*v+1]=0; } previous[(NUM_PREVIOUS+3)*v+2]=i; } void Ptngc_comp_to_lz77(unsigned int *vals, int nvals, - unsigned int *data, int *ndata, - unsigned int *len, int *nlens, - unsigned int *offsets, int *noffsets) + unsigned int *data, int *ndata, + unsigned int *len, int *nlens, + unsigned int *offsets, int *noffsets) { int noff=0; int ndat=0; @@ -192,112 +192,112 @@ void Ptngc_comp_to_lz77(unsigned int *vals, int nvals, #endif int firstoffset=i-MAX_OFFSET; if (firstoffset<0) - firstoffset=0; + firstoffset=0; if (i!=0) - { - int largest_len=0; - int largest_offset=0; - int icirc, ncirc; - /* Is this identical to a previous offset? Prefer close - values for offset. Search through circular buffer for the - possible values for the start of this string. */ - ncirc=previous[(NUM_PREVIOUS+3)*vals[i]]; - for (icirc=0; icirc=firstoffset) - { - for (k=0; i+klargest_len) && ((k>=(i-j)+16) || ((k>4) && (i-j==1)))) - { - largest_len=k; - largest_offset=j; - } - } - j++; - } - } + while ((j=firstoffset) + { + for (k=0; i+klargest_len) && ((k>=(i-j)+16) || ((k>4) && (i-j==1)))) + { + largest_len=k; + largest_offset=j; + } + } + j++; + } + } #if 0 - /* Search in sorted string buffer too. */ - kmin=info[i*2+1]-MAX_STRING_SEARCH; - kmax=info[i*2+1]+MAX_STRING_SEARCH; - if (kmin<0) - kmin=0; - if (kmax>=nvals) - kmax=nvals; - for (k=kmin; k=i)) - { - for (m=0; i+mlargest_len) && (m>4) && (m+2>=(i-s))) - { - largest_len=m; - largest_offset=s; + /* Search in sorted string buffer too. */ + kmin=info[i*2+1]-MAX_STRING_SEARCH; + kmax=info[i*2+1]+MAX_STRING_SEARCH; + if (kmin<0) + kmin=0; + if (kmax>=nvals) + kmax=nvals; + for (k=kmin; k=i)) + { + for (m=0; i+mlargest_len) && (m>4) && (m+2>=(i-s))) + { + largest_len=m; + largest_offset=s; #if 0 - fprintf(stderr,"Offset: %d %d\n",m,i-s); + fprintf(stderr,"Offset: %d %d\n",m,i-s); #endif - } - } - } + } + } + } #endif - /* Check how to write this info. */ - if (largest_len>MAX_LEN) - largest_len=MAX_LEN; - if (largest_len) - { - if (i-largest_offset==1) - { - data[ndat++]=0; - } - else - { - data[ndat++]=1; - offsets[noff++]=i-largest_offset; - } - len[nlen++]=largest_len; + /* Check how to write this info. */ + if (largest_len>MAX_LEN) + largest_len=MAX_LEN; + if (largest_len) + { + if (i-largest_offset==1) + { + data[ndat++]=0; + } + else + { + data[ndat++]=1; + offsets[noff++]=i-largest_offset; + } + len[nlen++]=largest_len; #if 0 - fprintf(stderr,"l:o: %d:%d data=%d i=%d\n",largest_len,i-largest_offset,ndat,i); - fflush(stderr); + fprintf(stderr,"l:o: %d:%d data=%d i=%d\n",largest_len,i-largest_offset,ndat,i); + fflush(stderr); #endif #if 0 - fprintf(stderr,"Found largest len %d at %d.\n",largest_len,i-largest_offset); + fprintf(stderr,"Found largest len %d at %d.\n",largest_len,i-largest_offset); #endif - /* Add these values to the circular buffer. */ - for (k=0; k=nvals) - { - fprintf(stderr,"too many vals.\n"); - exit(EXIT_FAILURE); - } - i++; - } - } + if (v==1) + offset=offsets[joff++]; + for (k=0; k=nvals) + { + fprintf(stderr,"too many vals.\n"); + exit(EXIT_FAILURE); + } + i++; + } + } else - vals[i++]=v-2; + vals[i++]=v-2; } } diff --git a/src/compression/merge_sort.c b/src/compression/merge_sort.c index 204243b..1d49903 100644 --- a/src/compression/merge_sort.c +++ b/src/compression/merge_sort.c @@ -18,12 +18,12 @@ #include "compression/merge_sort.h" static void ms_inner(void *base, size_t size, - size_t start, size_t end, - int (*compar)(const void *v1,const void *v2,const void *private), - const void *private, - char *workarray) + size_t start, size_t end, + int (*compar)(const void *v1,const void *v2,const void *private), + const void *private, + char *workarray) { - int middle; + size_t middle; if ((end-start)>1) { char *cbase=(char *)base; @@ -32,59 +32,59 @@ static void ms_inner(void *base, size_t size, printf("For start %d end %d obtained new middle: %d\n",start,end,middle); #endif ms_inner(base,size, - start,middle, - compar,private,workarray); + start,middle, + compar,private,workarray); ms_inner(base,size, - middle,end, - compar,private,workarray); + middle,end, + compar,private,workarray); #if 0 printf("For start %d end %d Before merge: Comparing element %d with %d\n",start,end,middle-1,middle); #endif if (compar(cbase+(middle-1)*size,cbase+middle*size,private)>0) - { - /* Merge to work array. */ - int i, n=end-start; - int ileft=start; - int iright=middle; - for (i=0; i0) - { - memcpy(workarray+i*size,cbase+iright*size,size); - iright++; - } - else - { - memcpy(workarray+i*size,cbase+ileft*size,size); - ileft++; - } - } - } - /* Copy result back. */ - memcpy(cbase+start*size,workarray,(end-start)*size); - } + { + /* Merge to work array. */ + size_t i, n=end-start; + size_t ileft=start; + size_t iright=middle; + for (i=0; i0) + { + memcpy(workarray+i*size,cbase+iright*size,size); + iright++; + } + else + { + memcpy(workarray+i*size,cbase+ileft*size,size); + ileft++; + } + } + } + /* Copy result back. */ + memcpy(cbase+start*size,workarray,(end-start)*size); + } } } void Ptngc_merge_sort(void *base, size_t nmemb, size_t size, - int (*compar)(const void *v1,const void *v2,const void *private), - void *private) + int (*compar)(const void *v1,const void *v2,const void *private), + void *private) { char *warr=warnmalloc(nmemb*size); ms_inner(base,size,0,nmemb,compar,private,warr); diff --git a/src/compression/mtf.c b/src/compression/mtf.c index a2ba3cf..1d74af3 100644 --- a/src/compression/mtf.c +++ b/src/compression/mtf.c @@ -19,7 +19,7 @@ /* Move to front coding. Acceptable inputs are max 8 bits (0-0xFF) */ static void comp_conv_to_mtf_byte(unsigned char *vals, int nvals, - unsigned char *valsmtf) + unsigned char *valsmtf) { int i; /* Indices into a linked list */ @@ -41,27 +41,27 @@ static void comp_conv_to_mtf_byte(unsigned char *vals, int nvals, int oldptr=-1; int r=0; while (dict[ptr]!=v) - { - oldptr=ptr; - ptr=list[ptr]; - r++; - } + { + oldptr=ptr; + ptr=list[ptr]; + r++; + } valsmtf[i]=(unsigned char)r; /* Move it to front in list */ /* Is it the head? Then it is already at the front. */ if (oldptr!=-1) - { - /* Remove it from inside the list */ - list[oldptr]=list[ptr]; - /* Move it to the front. */ - list[ptr]=head; - head=ptr; - } + { + /* Remove it from inside the list */ + list[oldptr]=list[ptr]; + /* Move it to the front. */ + list[ptr]=head; + head=ptr; + } } } void Ptngc_comp_conv_to_mtf_partial(unsigned int *vals, int nvals, - unsigned int *valsmtf) + unsigned int *valsmtf) { unsigned char *tmp=warnmalloc(nvals*2); int i, j; @@ -70,23 +70,23 @@ void Ptngc_comp_conv_to_mtf_partial(unsigned int *vals, int nvals, for (j=0; j<3; j++) { for (i=0; i>(8*j))&0xFF); + tmp[i]=(unsigned char)((vals[i]>>(8*j))&0xFF); comp_conv_to_mtf_byte(tmp,nvals,tmp+nvals); for (i=0; i>(8*j))&0xFF); + tmp[i]=(unsigned char)((vals[i]>>(8*j))&0xFF); comp_conv_to_mtf_byte(tmp,nvals,valsmtf+j*nvals); } free(tmp); @@ -94,7 +94,7 @@ void Ptngc_comp_conv_to_mtf_partial3(unsigned int *vals, int nvals, /* Move to front decoding */ static void comp_conv_from_mtf_byte(unsigned char *valsmtf, int nvals, - unsigned char *vals) + unsigned char *vals) { int i; /* Indices into a linked list */ @@ -116,27 +116,27 @@ static void comp_conv_from_mtf_byte(unsigned char *valsmtf, int nvals, int oldptr=-1; int cnt=0; while (cnt>(8*j))&0xFF); + tmp[i]=(unsigned char)((valsmtf[i]>>(8*j))&0xFF); comp_conv_from_mtf_byte(tmp,nvals,tmp+nvals); for (i=0; imin_rle) { /* Insert run-length */ unsigned int run=((unsigned int)nsim); while (run>1) - { - if (run&0x1U) - rle[(*j)++]=1; - else - rle[(*j)++]=0; - run>>=1; - } + { + if (run&0x1U) + rle[(*j)++]=1; + else + rle[(*j)++]=0; + run>>=1; + } nsim=1; } while (nsim--) @@ -39,8 +39,8 @@ static void add_rle(unsigned int *rle, Acceptable inputs are about 16 bits (0-0xFFFF) If input is 0-N output will be be values of 0-(N+2) */ void Ptngc_comp_conv_to_rle(unsigned int *vals, int nvals, - unsigned int *rle, int *nrle, - int min_rle) + unsigned int *rle, int *nrle, + int min_rle) { int i; int j=0; @@ -49,21 +49,21 @@ void Ptngc_comp_conv_to_rle(unsigned int *vals, int nvals, for (i=0; i=MAX_FVAL) - goto error; + if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL) + goto error; return 0; error: #if 0 for (iframe=0; iframe=MAX_FVAL) - 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_FVAL); + if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL) + 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_FVAL); #endif return 1; } @@ -56,46 +56,46 @@ static int verify_input_data_float(float *x, int natoms, int nframes, float prec for (iframe=0; iframe=MAX_FVAL) - goto error; + if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL) + goto error; return 0; error: #if 0 for (iframe=0; iframe=MAX_FVAL) - 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_FVAL); + if (fabs(x[iframe*natoms*3+i*3+j]/precision+0.5)>=MAX_FVAL) + 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_FVAL); #endif return 1; } static int quantize(double *x, int natoms, int nframes, - double precision, - int *quant) + double precision, + int *quant) { int iframe, i, j; for (iframe=0; iframe=2) - { - current_coding=TNG_COMPRESS_ALGO_POS_XTC3; - current_coding_parameter=0; - compress_quantized_pos(quant,NULL,quant_intra,natoms,1,speed, - current_coding,current_coding_parameter, - 0,0,prec_hi,prec_lo,¤t_code_size,NULL); - if (current_code_size=6) - { - current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTRA; - current_coding_parameter=0; - compress_quantized_pos(quant,NULL,quant_intra,natoms,1,speed, - current_coding,current_coding_parameter, - 0,0,prec_hi,prec_lo,¤t_code_size,NULL); - if (current_code_size=4) - { - current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTER; - current_coding_parameter=0; - compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed, - TNG_COMPRESS_ALGO_POS_XTC2,0, - current_coding,current_coding_parameter, - prec_hi,prec_lo,¤t_code_size,NULL); - current_code_size-=initial_code_size; /* Correct for the use of XTC2 for the first frame. */ - if (current_code_size=6) - { - current_coding=TNG_COMPRESS_ALGO_POS_BWLZH_INTRA; - current_coding_parameter=0; - compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed, - TNG_COMPRESS_ALGO_POS_XTC2,0, - current_coding,current_coding_parameter, - prec_hi,prec_lo,¤t_code_size,NULL); - current_code_size-=initial_code_size; /* Correct for the use of XTC2 for the first frame. */ - if (current_code_size=4) - { - current_coding=TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE; - current_coding_parameter=0; - compress_quantized_vel(quant,NULL,natoms,1,speed, - current_coding,current_coding_parameter, - 0,0,prec_hi,prec_lo,¤t_code_size,NULL); - if ((best_coding==-1) || (current_code_size=4) - { - /* Test BWLZH inter */ - current_coding=TNG_COMPRESS_ALGO_VEL_BWLZH_INTER; - current_coding_parameter=0; - compress_quantized_vel(quant,quant_inter,natoms,nframes,speed, - TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE,initial_numbits, - current_coding,current_coding_parameter, - prec_hi,prec_lo,¤t_code_size,NULL); - current_code_size-=initial_code_size; /* Correct for the initial frame */ - if (current_code_size1) { if (coding==-1) - { - coding_parameter=-1; - determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo, - &coding,&coding_parameter); - } + { + coding_parameter=-1; + determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo, + &coding,&coding_parameter); + } else if (coding_parameter==-1) - { - determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo, - &coding,&coding_parameter); - } + { + determine_best_pos_coding(quant,quant_inter,quant_intra,natoms,nframes,speed,prec_hi,prec_lo, + &coding,&coding_parameter); + } } compress_quantized_pos(quant,quant_inter,quant_intra,natoms,nframes,speed, - initial_coding,initial_coding_parameter, - coding,coding_parameter, - prec_hi,prec_lo,nitems,data); + initial_coding,initial_coding_parameter, + coding,coding_parameter, + prec_hi,prec_lo,nitems,data); free(quant_inter); free(quant_intra); if (algo[0]==-1) @@ -1217,9 +1217,9 @@ char DECLSPECDLLEXPORT *tng_compress_pos_int(int *pos, int natoms, int nframes, } char DECLSPECDLLEXPORT *tng_compress_pos(double *pos, int natoms, int nframes, - double desired_precision, - int speed,int *algo, - int *nitems) + double desired_precision, + int speed,int *algo, + int *nitems) { int *quant=malloc(natoms*nframes*3*sizeof *quant); char *data; @@ -1235,9 +1235,9 @@ char DECLSPECDLLEXPORT *tng_compress_pos(double *pos, int natoms, int nframes, } char DECLSPECDLLEXPORT *tng_compress_pos_float(float *pos, int natoms, int nframes, - float desired_precision, - int speed,int *algo, - int *nitems) + float desired_precision, + int speed,int *algo, + int *nitems) { int *quant=malloc(natoms*nframes*3*sizeof *quant); char *data; @@ -1253,10 +1253,10 @@ char DECLSPECDLLEXPORT *tng_compress_pos_float(float *pos, int natoms, int nfram } char DECLSPECDLLEXPORT *tng_compress_pos_find_algo(double *pos, int natoms, int nframes, - double desired_precision, - int speed, - int *algo, - int *nitems) + double desired_precision, + int speed, + int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1266,10 +1266,10 @@ char DECLSPECDLLEXPORT *tng_compress_pos_find_algo(double *pos, int natoms, int } char DECLSPECDLLEXPORT *tng_compress_pos_float_find_algo(float *pos, int natoms, int nframes, - float desired_precision, - int speed, - int *algo, - int *nitems) + float desired_precision, + int speed, + int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1279,9 +1279,9 @@ char DECLSPECDLLEXPORT *tng_compress_pos_float_find_algo(float *pos, int natoms, } char DECLSPECDLLEXPORT *tng_compress_pos_int_find_algo(int *pos, int natoms, int nframes, - unsigned long prec_hi, unsigned long prec_lo, - int speed,int *algo, - int *nitems) + unsigned long prec_hi, unsigned long prec_lo, + int speed,int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1303,13 +1303,13 @@ int DECLSPECDLLEXPORT tng_compress_nalgo(void) } char DECLSPECDLLEXPORT *tng_compress_vel_int(int *vel, int natoms, int nframes, - unsigned long prec_hi, unsigned long prec_lo, - int speed, int *algo, - int *nitems) + unsigned long prec_hi, unsigned long prec_lo, + int speed, int *algo, + int *nitems) { char *data=malloc(natoms*nframes*14+11*4); /* 12 bytes are required to store 4 32 bit integers - This is 17% extra. The final 11*4 is to store information - needed for decompression. */ + This is 17% extra. The final 11*4 is to store information + needed for decompression. */ int *quant=vel; int *quant_inter=malloc(natoms*nframes*3*sizeof *quant_inter); @@ -1334,12 +1334,12 @@ char DECLSPECDLLEXPORT *tng_compress_vel_int(int *vel, int natoms, int nframes, { initial_coding_parameter=-1; determine_best_vel_initial_coding(quant,natoms,speed,prec_hi,prec_lo, - &initial_coding,&initial_coding_parameter); + &initial_coding,&initial_coding_parameter); } else if (initial_coding_parameter==-1) { determine_best_vel_initial_coding(quant,natoms,speed,prec_hi,prec_lo, - &initial_coding,&initial_coding_parameter); + &initial_coding,&initial_coding_parameter); } if (nframes==1) @@ -1351,22 +1351,22 @@ char DECLSPECDLLEXPORT *tng_compress_vel_int(int *vel, int natoms, int nframes, if (nframes>1) { if (coding==-1) - { - coding_parameter=-1; - determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo, - &coding,&coding_parameter); - } + { + coding_parameter=-1; + determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo, + &coding,&coding_parameter); + } else if (coding_parameter==-1) - { - determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo, - &coding,&coding_parameter); - } + { + determine_best_vel_coding(quant,quant_inter,natoms,nframes,speed,prec_hi,prec_lo, + &coding,&coding_parameter); + } } compress_quantized_vel(quant,quant_inter,natoms,nframes,speed, - initial_coding,initial_coding_parameter, - coding,coding_parameter, - prec_hi,prec_lo,nitems,data); + initial_coding,initial_coding_parameter, + coding,coding_parameter, + prec_hi,prec_lo,nitems,data); free(quant_inter); if (algo[0]==-1) algo[0]=initial_coding; @@ -1380,9 +1380,9 @@ char DECLSPECDLLEXPORT *tng_compress_vel_int(int *vel, int natoms, int nframes, } char DECLSPECDLLEXPORT *tng_compress_vel(double *vel, int natoms, int nframes, - double desired_precision, - int speed, int *algo, - int *nitems) + double desired_precision, + int speed, int *algo, + int *nitems) { int *quant=malloc(natoms*nframes*3*sizeof *quant); char *data; @@ -1397,9 +1397,9 @@ char DECLSPECDLLEXPORT *tng_compress_vel(double *vel, int natoms, int nframes, } char DECLSPECDLLEXPORT *tng_compress_vel_float(float *vel, int natoms, int nframes, - float desired_precision, - int speed, int *algo, - int *nitems) + float desired_precision, + int speed, int *algo, + int *nitems) { int *quant=malloc(natoms*nframes*3*sizeof *quant); char *data; @@ -1414,10 +1414,10 @@ char DECLSPECDLLEXPORT *tng_compress_vel_float(float *vel, int natoms, int nfram } char DECLSPECDLLEXPORT *tng_compress_vel_find_algo(double *vel, int natoms, int nframes, - double desired_precision, - int speed, - int *algo, - int *nitems) + double desired_precision, + int speed, + int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1427,10 +1427,10 @@ char DECLSPECDLLEXPORT *tng_compress_vel_find_algo(double *vel, int natoms, int } char DECLSPECDLLEXPORT *tng_compress_vel_float_find_algo(float *vel, int natoms, int nframes, - float desired_precision, - int speed, - int *algo, - int *nitems) + float desired_precision, + int speed, + int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1440,10 +1440,10 @@ char DECLSPECDLLEXPORT *tng_compress_vel_float_find_algo(float *vel, int natoms, } char DECLSPECDLLEXPORT *tng_compress_vel_int_find_algo(int *vel, int natoms, int nframes, - unsigned long prec_hi, unsigned long prec_lo, - int speed, - int *algo, - int *nitems) + unsigned long prec_hi, unsigned long prec_lo, + int speed, + int *algo, + int *nitems) { algo[0]=-1; algo[1]=-1; @@ -1453,8 +1453,8 @@ char DECLSPECDLLEXPORT *tng_compress_vel_int_find_algo(int *vel, int natoms, int } int DECLSPECDLLEXPORT tng_compress_inquire(char *data,int *vel, int *natoms, - int *nframes, double *precision, - int *algo) + int *nframes, double *precision, + int *algo) { int bufloc=0; fix_t prec_hi, prec_lo; @@ -1549,7 +1549,7 @@ static int tng_compress_uncompress_pos_gen(char *data,double *posd,float *posf,i /* The initial frame */ coder=Ptngc_coder_init(); rval=Ptngc_unpack_array(coder,(unsigned char*)data+bufloc,quant,natoms*3, - initial_coding,initial_coding_parameter,natoms); + initial_coding,initial_coding_parameter,natoms); Ptngc_coder_deinit(coder); if (rval) goto error; @@ -1561,21 +1561,21 @@ static int tng_compress_uncompress_pos_gen(char *data,double *posd,float *posf,i (initial_coding==TNG_COMPRESS_ALGO_POS_XTC3)) { if (posd) - unquantize(posd,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); + unquantize(posd,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); else if (posf) - unquantize_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); + unquantize_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); else if (posi) - memcpy(posi,quant,natoms*3*sizeof *posi); + memcpy(posi,quant,natoms*3*sizeof *posi); } else if ((initial_coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA) || - (initial_coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA)) + (initial_coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA)) { if (posd) - unquantize_intra_differences(posd,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); + unquantize_intra_differences(posd,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); else if (posf) - unquantize_intra_differences_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); + unquantize_intra_differences_float(posf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); else if (posi) - unquantize_intra_differences_int(posi,natoms,1,quant); + unquantize_intra_differences_int(posi,natoms,1,quant); unquant_intra_differences_first_frame(quant,natoms); } /* The remaining frames. */ @@ -1584,45 +1584,45 @@ static int tng_compress_uncompress_pos_gen(char *data,double *posd,float *posf,i bufloc+=4; coder=Ptngc_coder_init(); rval=Ptngc_unpack_array(coder,(unsigned char *)data+bufloc,quant+natoms*3,(nframes-1)*natoms*3, - coding,coding_parameter,natoms); + coding,coding_parameter,natoms); Ptngc_coder_deinit(coder); if (rval) - goto error; + goto error; if ((coding==TNG_COMPRESS_ALGO_POS_STOPBIT_INTER) || - (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTER) || - (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTER)) - { - /* This requires that the first frame is already in one-to-one format, even if intra-frame - compression was done there. Therefore the unquant_intra_differences_first_frame should be called - before to convert it correctly. */ - if (posd) - unquantize_inter_differences(posd,natoms,nframes,PRECISION(*prec_hi,*prec_lo),quant); - else if (posf) - unquantize_inter_differences_float(posf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo),quant); - else if (posi) - unquantize_inter_differences_int(posi,natoms,nframes,quant); - } + (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTER) || + (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTER)) + { + /* This requires that the first frame is already in one-to-one format, even if intra-frame + compression was done there. Therefore the unquant_intra_differences_first_frame should be called + before to convert it correctly. */ + if (posd) + unquantize_inter_differences(posd,natoms,nframes,PRECISION(*prec_hi,*prec_lo),quant); + else if (posf) + unquantize_inter_differences_float(posf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo),quant); + else if (posi) + unquantize_inter_differences_int(posi,natoms,nframes,quant); + } else if ((coding==TNG_COMPRESS_ALGO_POS_XTC2) || - (coding==TNG_COMPRESS_ALGO_POS_XTC3) || - (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE)) - { - if (posd) - unquantize(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (posf) - unquantize_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (posi) - memcpy(posi+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *posi); - } + (coding==TNG_COMPRESS_ALGO_POS_XTC3) || + (coding==TNG_COMPRESS_ALGO_POS_TRIPLET_ONETOONE)) + { + if (posd) + unquantize(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (posf) + unquantize_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (posi) + memcpy(posi+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *posi); + } else if ((coding==TNG_COMPRESS_ALGO_POS_TRIPLET_INTRA) || - (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA)) - { - if (posd) - unquantize_intra_differences(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (posf) - unquantize_intra_differences_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (posi) - unquantize_intra_differences_int(posi+natoms*3,natoms,nframes-1,quant+natoms*3); - } + (coding==TNG_COMPRESS_ALGO_POS_BWLZH_INTRA)) + { + if (posd) + unquantize_intra_differences(posd+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (posf) + unquantize_intra_differences_float(posf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (posi) + unquantize_intra_differences_int(posi+natoms*3,natoms,nframes-1,quant+natoms*3); + } } error: free(quant); @@ -1696,7 +1696,7 @@ static int tng_compress_uncompress_vel_gen(char *data,double *veld,float *velf,i /* The initial frame */ coder=Ptngc_coder_init(); rval=Ptngc_unpack_array(coder,(unsigned char*)data+bufloc,quant,natoms*3, - initial_coding,initial_coding_parameter,natoms); + initial_coding,initial_coding_parameter,natoms); Ptngc_coder_deinit(coder); if (rval) goto error; @@ -1708,11 +1708,11 @@ static int tng_compress_uncompress_vel_gen(char *data,double *veld,float *velf,i (initial_coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE)) { if (veld) - unquantize(veld,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); + unquantize(veld,natoms,1,PRECISION(*prec_hi,*prec_lo),quant); else if (velf) - unquantize_float(velf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); + unquantize_float(velf,natoms,1,(float)PRECISION(*prec_hi,*prec_lo),quant); else if (veli) - memcpy(veli,quant,natoms*3*sizeof *veli); + memcpy(veli,quant,natoms*3*sizeof *veli); } /* The remaining frames. */ if (nframes>1) @@ -1720,35 +1720,35 @@ static int tng_compress_uncompress_vel_gen(char *data,double *veld,float *velf,i bufloc+=4; coder=Ptngc_coder_init(); rval=Ptngc_unpack_array(coder,(unsigned char *)data+bufloc,quant+natoms*3,(nframes-1)*natoms*3, - coding,coding_parameter,natoms); + coding,coding_parameter,natoms); Ptngc_coder_deinit(coder); if (rval) - goto error; + goto error; /* Inter-frame compression? */ if ((coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_INTER) || - (coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER) || - (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_INTER)) - { - /* This requires that the first frame is already in one-to-one format. */ - if (veld) - unquantize_inter_differences(veld,natoms,nframes,PRECISION(*prec_hi,*prec_lo),quant); - else if (velf) - unquantize_inter_differences_float(velf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo),quant); - else if (veli) - unquantize_inter_differences_int(veli,natoms,nframes,quant); - } + (coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_INTER) || + (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_INTER)) + { + /* This requires that the first frame is already in one-to-one format. */ + if (veld) + unquantize_inter_differences(veld,natoms,nframes,PRECISION(*prec_hi,*prec_lo),quant); + else if (velf) + unquantize_inter_differences_float(velf,natoms,nframes,(float)PRECISION(*prec_hi,*prec_lo),quant); + else if (veli) + unquantize_inter_differences_int(veli,natoms,nframes,quant); + } /* One-to-one compression? */ else if ((coding==TNG_COMPRESS_ALGO_VEL_STOPBIT_ONETOONE) || - (coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE) || - (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE)) - { - if (veld) - unquantize(veld+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (velf) - unquantize_float(velf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); - else if (veli) - memcpy(veli+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *veli); - } + (coding==TNG_COMPRESS_ALGO_VEL_TRIPLET_ONETOONE) || + (coding==TNG_COMPRESS_ALGO_VEL_BWLZH_ONETOONE)) + { + if (veld) + unquantize(veld+natoms*3,natoms,nframes-1,PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (velf) + unquantize_float(velf+natoms*3,natoms,nframes-1,(float)PRECISION(*prec_hi,*prec_lo),quant+natoms*3); + else if (veli) + memcpy(veli+natoms*3,quant+natoms*3,natoms*3*(nframes-1)*sizeof *veli); + } } error: free(quant); @@ -1811,15 +1811,15 @@ int DECLSPECDLLEXPORT tng_compress_uncompress_int(char *data,int *posvel, unsign } void DECLSPECDLLEXPORT tng_compress_int_to_double(int *posvel_int,unsigned long prec_hi, unsigned long prec_lo, - int natoms,int nframes, - double *posvel_double) + int natoms,int nframes, + double *posvel_double) { unquantize(posvel_double,natoms,nframes,PRECISION(prec_hi,prec_lo),posvel_int); } void DECLSPECDLLEXPORT tng_compress_int_to_float(int *posvel_int,unsigned long prec_hi, unsigned long prec_lo, - int natoms,int nframes, - float *posvel_float) + int natoms,int nframes, + float *posvel_float) { unquantize_float(posvel_float,natoms,nframes,(float)PRECISION(prec_hi,prec_lo),posvel_int); } diff --git a/src/compression/vals16.c b/src/compression/vals16.c index df44d8d..10a1cd8 100644 --- a/src/compression/vals16.c +++ b/src/compression/vals16.c @@ -16,29 +16,29 @@ /* Coding 32 bit ints in sequences of 16 bit ints. Worst case the output is 3*nvals long. */ void Ptngc_comp_conv_to_vals16(unsigned int *vals,int nvals, - unsigned int *vals16, int *nvals16) + unsigned int *vals16, int *nvals16) { int i; int j=0; for (i=0; i>15; - vals16[j++]=lo; - if (hi<=0x7FFFU) - vals16[j++]=hi; - else - { - unsigned int lohi=(hi&0x7FFFU)|0x8000U; - unsigned int hihi=hi>>15; - vals16[j++]=lohi; - vals16[j++]=hihi; - } - } + { + unsigned int lo=(vals[i]&0x7FFFU)|0x8000U; + unsigned int hi=vals[i]>>15; + vals16[j++]=lo; + if (hi<=0x7FFFU) + vals16[j++]=hi; + else + { + unsigned int lohi=(hi&0x7FFFU)|0x8000U; + unsigned int hihi=hi>>15; + vals16[j++]=lohi; + vals16[j++]=hihi; + } + } } #if 0 /* Test that things that detect that this is bad really works. */ @@ -48,26 +48,26 @@ void Ptngc_comp_conv_to_vals16(unsigned int *vals,int nvals, } void Ptngc_comp_conv_from_vals16(unsigned int *vals16,int nvals16, - unsigned int *vals, int *nvals) + unsigned int *vals, int *nvals) { int i=0; int j=0; while (i>bits; - + x=a_U*b_L; x_LU+=x & L_m; x_UL=x>>bits; - + x=a_L*b_U; x_LU+=x & L_m; x_UL+=x>>bits; @@ -105,9 +105,9 @@ void Ptngc_widediv(unsigned int hi, unsigned int lo, unsigned int i, unsigned in { #if defined(TRAJNG_X86_GCC_INLINE_MULDIV) __asm__ __volatile__ ("divl %%ecx\n\t" - : "=a" (lo), "=d" (hi) - : "a" (lo),"d" (hi), "c" (i) - : "cc"); + : "=a" (lo), "=d" (hi) + : "a" (lo),"d" (hi), "c" (i) + : "cc"); *result=lo; *remainder=hi; #else /* TRAJNG X86 GCC INLINE MULDIV */ @@ -126,7 +126,7 @@ void Ptngc_widediv(unsigned int hi, unsigned int lo, unsigned int i, unsigned in unsigned int hibit=bits2-1U; unsigned int hibit_mask=1U<lo) - { - unsigned int t; - hi--; /* Borrow */ - t=allbits-s_L; - lo+=t+1; - } - else - lo-=s_L; - - /* Set bit. */ - res|=rmask; - } + { + /* Subtract */ + hi-=s_U; + if (s_L>lo) + { + unsigned int t; + hi--; /* Borrow */ + t=allbits-s_L; + lo+=t+1; + } + else + lo-=s_L; + + /* Set bit. */ + res|=rmask; + } rmask>>=1; s_L>>=1; if (s_U & 1U) - s_L|=hibit_mask; + s_L|=hibit_mask; s_U>>=1; } *remainder=lo; @@ -179,7 +179,7 @@ static void largeint_add_gen(unsigned int v1, unsigned int *largeint, int n, int v2=(largeint[j]+carry)&0xFFFFFFFFU; carry=0; if ((((unsigned int)-1)&0xFFFFFFFFU) -164 mul */ - largeint_add_gen(lo,largeint_out,n,i); - if (i+164 mul */ + largeint_add_gen(lo,largeint_out,n,i); + if (i+1>=1; if (!extract_mask) - { - extract_mask=0x80U; - *ptr=(*ptr)+1; - *bitptr=0; - if (nbits) - thisval=**ptr; - } + { + extract_mask=0x80U; + *ptr=(*ptr)+1; + *bitptr=0; + if (nbits) + thisval=**ptr; + } } #ifdef SHOWIT fprintf(stderr,"%x\n",val); @@ -301,25 +301,25 @@ static int read_instruction(unsigned char **ptr, int *bitptr) { bits=readbits(ptr,bitptr,1); if (!bits) - instr=INSTR_BASE_RUNLENGTH; + instr=INSTR_BASE_RUNLENGTH; else - { - bits=readbits(ptr,bitptr,2); - if (bits==0) - instr=INSTR_ONLY_LARGE; - else if (bits==1) - instr=INSTR_ONLY_SMALL; - else if (bits==2) - instr=INSTR_LARGE_BASE_CHANGE; - else if (bits==3) - { - bits=readbits(ptr,bitptr,1); - if (bits==0) - instr=INSTR_FLIP; - else - instr=INSTR_LARGE_RLE; - } - } + { + bits=readbits(ptr,bitptr,2); + if (bits==0) + instr=INSTR_ONLY_LARGE; + else if (bits==1) + instr=INSTR_ONLY_SMALL; + else if (bits==2) + instr=INSTR_LARGE_BASE_CHANGE; + else if (bits==3) + { + bits=readbits(ptr,bitptr,1); + if (bits==0) + instr=INSTR_FLIP; + else + instr=INSTR_LARGE_RLE; + } + } } return instr; } @@ -346,12 +346,12 @@ static void swap_is_better(int *input, int *minint, int *sum_normal, int *sum_sw normal[2]=input[6+i]-input[3+i]; /* minint[i]-minint[i] cancels out */ swap_ints(normal,swapped); for (j=1; j<3; j++) - { - if (positive_int(normal[j])>(unsigned int)normal_max) - normal_max=positive_int(normal[j]); - if (positive_int(swapped[j])>(unsigned int)swapped_max) - swapped_max=positive_int(swapped[j]); - } + { + if (positive_int(normal[j])>(unsigned int)normal_max) + normal_max=positive_int(normal[j]); + if (positive_int(swapped[j])>(unsigned int)swapped_max) + swapped_max=positive_int(swapped[j]); + } } if (normal_max==0) normal_max=1; @@ -381,21 +381,21 @@ static void swapdecide(struct coder *coder, int *input,int *swapatoms, int *larg ((normal>shift) & 0xFFU; - shift+=8; - } + { + result[i*4+j]=(largeint[i]>>shift) & 0xFFU; + shift+=8; + } } } @@ -494,10 +494,10 @@ static void trajcoder_base_decompress(unsigned char *input, int n, int *index, i int shift=0; largeint[i]=0U; for (j=0; j<4; j++) - { - largeint[i]|=((unsigned int)input[i*4+j])<magic[small_index+QUITE_LARGE]) - { - is=1; - break; - } + if (positive_int(input[i])>magic[small_index+QUITE_LARGE]) + { + is=1; + break; + } } return is; } @@ -574,21 +574,21 @@ static void insert_batch(int *input_ptr, int ntriplets_left, int *prevcoord,int { int i; for (i=0; imaxint[j]) - maxint[j]=input[i*3+j]; - if (input[i*3+j]maxint[j]) + maxint[j]=input[i*3+j]; + if (input[i*3+j]intmax) - if (slargest_required_base) - largest_required_base=encode_ints[ienc]; - /* Also compute what the largest base is for the current runlength setting! */ - largest_runlength_base=0; - for (ienc=0; (ienclargest_runlength_base) - largest_runlength_base=encode_ints[ienc]; - - largest_required_index=Ptngc_find_magic_index(largest_required_base); - largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); - - if (largest_required_indexntriplets_left) - new_runlength=ntriplets_left; - - /* We must at least try to get some small integers going. */ - if (new_runlength==0) - { - new_runlength=1; - new_small_index=small_index; - } - - iter_runlength=new_runlength; - iter_small_index=new_small_index; - - /* Iterate to find optimal encoding and runlength */ + fprintf(stderr,"Update batch due to large int.\n"); +#endif + if ((swapatoms) && (didswap)) + { + /* Keep swapped values. */ + for (i=0; i<2; i++) + for (ienc=0; ienc<3; ienc++) + encode_ints[i*3+ienc]=encode_ints[(i+1)*3+ienc]; + } + insert_batch(input_ptr,ntriplets_left,prevcoord,minint,encode_ints,min_runlength,&nencode); + } + /* Here we should only have differences for the atom coordinates. */ + /* Convert the ints to positive ints */ + for (ienc=0; ienclargest_required_base) + largest_required_base=encode_ints[ienc]; + /* Also compute what the largest base is for the current runlength setting! */ + largest_runlength_base=0; + for (ienc=0; (ienclargest_runlength_base) + largest_runlength_base=encode_ints[ienc]; + + largest_required_index=Ptngc_find_magic_index(largest_required_base); + largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); + + if (largest_required_indexntriplets_left) + new_runlength=ntriplets_left; + + /* We must at least try to get some small integers going. */ + if (new_runlength==0) + { + new_runlength=1; + new_small_index=small_index; + } + + iter_runlength=new_runlength; + iter_small_index=new_small_index; + + /* Iterate to find optimal encoding and runlength */ #ifdef SHOWIT - fprintf(stderr,"Entering iterative loop.\n"); - fflush(stderr); + fprintf(stderr,"Entering iterative loop.\n"); + fflush(stderr); #endif - do { - new_runlength=iter_runlength; - new_small_index=iter_small_index; + do { + new_runlength=iter_runlength; + new_small_index=iter_small_index; #ifdef SHOWIT - fprintf(stderr,"Test new_small_index=%d Base=%d\n",new_small_index,magic[new_small_index]); -#endif - /* What is the largest runlength - we can do with the currently - selected encoding? Also the max supported runlength is 6 triplets! */ - for (ienc=0; iencnew_small_index) - break; - } - if (ienc/3>new_runlength) - { - iter_runlength=ienc/3; + fprintf(stderr,"Test new_small_index=%d Base=%d\n",new_small_index,magic[new_small_index]); +#endif + /* What is the largest runlength + we can do with the currently + selected encoding? Also the max supported runlength is 6 triplets! */ + for (ienc=0; iencnew_small_index) + break; + } + if (ienc/3>new_runlength) + { + iter_runlength=ienc/3; #ifdef SHOWIT - fprintf(stderr,"I found a new possible runlength: %d\n",iter_runlength); -#endif - } - - /* How large encoding do we have to use? */ - largest_runlength_base=0; - for (ienc=0; ienclargest_runlength_base) - largest_runlength_base=encode_ints[ienc]; - largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); - if (largest_runlength_index!=new_small_index) - { - iter_small_index=largest_runlength_index; + fprintf(stderr,"I found a new possible runlength: %d\n",iter_runlength); +#endif + } + + /* How large encoding do we have to use? */ + largest_runlength_base=0; + for (ienc=0; ienclargest_runlength_base) + largest_runlength_base=encode_ints[ienc]; + largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); + if (largest_runlength_index!=new_small_index) + { + iter_small_index=largest_runlength_index; #ifdef SHOWIT - fprintf(stderr,"I found a new possible small index: %d Base=%d\n",iter_small_index,magic[iter_small_index]); + fprintf(stderr,"I found a new possible small index: %d Base=%d\n",iter_small_index,magic[iter_small_index]); #endif - } - } while ((new_runlength!=iter_runlength) || - (new_small_index!=iter_small_index)); + } + } while ((new_runlength!=iter_runlength) || + (new_small_index!=iter_small_index)); #ifdef SHOWIT - fprintf(stderr,"Exit iterative loop.\n"); - fflush(stderr); -#endif - - /* Verify that we got something good. We may have caught a - substantially larger atom. If so we should just bail - out and let the loop get on another lap. We may have a - minimum runlength though and then we have to fulfill - the request to write out these atoms! */ - rle_index_dep=0; - if (new_runlength<3) - rle_index_dep=IS_LARGE; - else if (new_runlength<6) - rle_index_dep=QUITE_LARGE; - if ((min_runlength) - || ((new_small_index=%g?\n",change,isum,(double)magic[small_index+change]*(double)magic[small_index+change]); + fprintf(stderr,"Tested decrease %d of index: %g>=%g?\n",change,isum,(double)magic[small_index+change]*(double)magic[small_index+change]); #endif - if (isum>(double)magic[small_index+change]*(double)magic[small_index+change]) - { + if (isum>(double)magic[small_index+change]*(double)magic[small_index+change]) + { #ifdef SHOWIT - fprintf(stderr,"Rejected decrease %d of index due to length of vector: %g>=%g\n",change,isum,(double)magic[small_index+change]*(double)magic[small_index+change]); -#endif - rejected=1; - change++; - } - } while ((change<0) && (rejected)); - if (change==0) - break; - } - } - - /* If the only thing to do is to change the base by - only one -1 it is probably not worth it. */ - if (!((change==-1) && (runlength==new_runlength))) - { - /* If we have a very short runlength we do not - want to do large base changes. It costs 6 - extra bits to do -2. We gain 2/3 - bits per value to decrease the index by -2, - ie 2 bits, so to any changes down we must - have a runlength of 3 or more to do it for - one molecule! If we have several molecules we - will gain of course, so not be so strict. */ - if ((change==-2) && (new_runlength<3)) - { - if (runlength==new_runlength) - change=0; - else - change=-1; + fprintf(stderr,"Rejected decrease %d of index due to length of vector: %g>=%g\n",change,isum,(double)magic[small_index+change]*(double)magic[small_index+change]); +#endif + rejected=1; + change++; + } + } while ((change<0) && (rejected)); + if (change==0) + break; + } + } + + /* If the only thing to do is to change the base by + only one -1 it is probably not worth it. */ + if (!((change==-1) && (runlength==new_runlength))) + { + /* If we have a very short runlength we do not + want to do large base changes. It costs 6 + extra bits to do -2. We gain 2/3 + bits per value to decrease the index by -2, + ie 2 bits, so to any changes down we must + have a runlength of 3 or more to do it for + one molecule! If we have several molecules we + will gain of course, so not be so strict. */ + if ((change==-2) && (new_runlength<3)) + { + if (runlength==new_runlength) + change=0; + else + change=-1; #ifdef SHOWIT - fprintf(stderr,"Rejected change by -2 due to too short runlenght. Change set to %d\n",change); -#endif - } - - /* First adjust base using large base change instruction (if necessary) */ - while ((change>1) || (change<-1) || ((new_runlength==6) && (change))) - { - unsigned int code=0U; - int this_change=change; - if (this_change>2) - this_change=2; - if (this_change<-2) - this_change=-2; - change-=this_change; + fprintf(stderr,"Rejected change by -2 due to too short runlenght. Change set to %d\n",change); +#endif + } + + /* First adjust base using large base change instruction (if necessary) */ + while ((change>1) || (change<-1) || ((new_runlength==6) && (change))) + { + unsigned int code=0U; + int this_change=change; + if (this_change>2) + this_change=2; + if (this_change<-2) + this_change=-2; + change-=this_change; #ifdef SHOWIT - fprintf(stderr,"Large base change: %d.\n",this_change); -#endif - small_index+=this_change; - if (this_change<0) - { - code|=2U; - this_change=-this_change; - } - code|=(unsigned int)(this_change-1); - write_instruction(coder,INSTR_LARGE_BASE_CHANGE,&output_ptr); - Ptngc_writebits(coder,code,2,&output_ptr); - } - /* If there is still base change or runlength changes to do, we do them now. */ - if ((new_runlength!=runlength) || (change)) - { - unsigned int ichange=(unsigned int)(change+1); - unsigned int code=0U; - unsigned int irun=(unsigned int)((new_runlength-1)*3); - if (new_runlength==6) - ichange=0; /* Means no change. The change has been taken care of explicitly by a large - base change instruction above. */ - code=ichange+irun; + fprintf(stderr,"Large base change: %d.\n",this_change); +#endif + small_index+=this_change; + if (this_change<0) + { + code|=2U; + this_change=-this_change; + } + code|=(unsigned int)(this_change-1); + write_instruction(coder,INSTR_LARGE_BASE_CHANGE,&output_ptr); + Ptngc_writebits(coder,code,2,&output_ptr); + } + /* If there is still base change or runlength changes to do, we do them now. */ + if ((new_runlength!=runlength) || (change)) + { + unsigned int ichange=(unsigned int)(change+1); + unsigned int code=0U; + unsigned int irun=(unsigned int)((new_runlength-1)*3); + if (new_runlength==6) + ichange=0; /* Means no change. The change has been taken care of explicitly by a large + base change instruction above. */ + code=ichange+irun; #ifdef SHOWIT - fprintf(stderr,"Small base change: %d Runlength change: %d\n",change,new_runlength); -#endif - small_index+=change; - write_instruction(coder,INSTR_BASE_RUNLENGTH,&output_ptr); - Ptngc_writebits(coder,code,4,&output_ptr); - runlength=new_runlength; - } - } + fprintf(stderr,"Small base change: %d Runlength change: %d\n",change,new_runlength); +#endif + small_index+=change; + write_instruction(coder,INSTR_BASE_RUNLENGTH,&output_ptr); + Ptngc_writebits(coder,code,4,&output_ptr); + runlength=new_runlength; + } + } #ifdef SHOWIT - else - fprintf(stderr,"Rejected base change due to only change==-1\n"); + else + fprintf(stderr,"Rejected base change due to only change==-1\n"); #endif #ifdef SHOWIT - fprintf(stderr,"Current small index: %d Base=%d\n",small_index,magic[small_index]); -#endif - } - /* If we have a large previous integer we can combine it with a sequence of small ints. */ - if (has_large) - { - /* If swapatoms is set to 1 but we did actually not - do any swapping, we must first write out the - large atom and then the small. If swapatoms is 1 - and we did swapping we can use the efficient - encoding. */ - if ((swapatoms) && (!didswap)) - { + fprintf(stderr,"Current small index: %d Base=%d\n",small_index,magic[small_index]); +#endif + } + /* If we have a large previous integer we can combine it with a sequence of small ints. */ + if (has_large) + { + /* If swapatoms is set to 1 but we did actually not + do any swapping, we must first write out the + large atom and then the small. If swapatoms is 1 + and we did swapping we can use the efficient + encoding. */ + if ((swapatoms) && (!didswap)) + { #ifdef SHOWIT - fprintf(stderr,"Swapatoms was set to 1 but we did not do swapping!\n"); - fprintf(stderr,"Only one large integer.\n"); + fprintf(stderr,"Swapatoms was set to 1 but we did not do swapping!\n"); + fprintf(stderr,"Only one large integer.\n"); #endif - /* Flush all large atoms. */ - flush_large(coder,&has_large,has_large_ints,has_large,large_index,large_nbits,compress_buffer,&output_ptr); + /* Flush all large atoms. */ + flush_large(coder,&has_large,has_large_ints,has_large,large_index,large_nbits,compress_buffer,&output_ptr); #ifdef SHOWIT - fprintf(stderr,"Sequence of only small integers.\n"); + fprintf(stderr,"Sequence of only small integers.\n"); #endif - write_instruction(coder,INSTR_ONLY_SMALL,&output_ptr); - } - else - { + write_instruction(coder,INSTR_ONLY_SMALL,&output_ptr); + } + else + { #ifdef SHOWIT - fprintf(stderr,"Sequence of one large and small integers (good compression).\n"); -#endif - /* Flush all large atoms but one! */ - if (has_large>1) - flush_large(coder,&has_large,has_large_ints,has_large-1,large_index,large_nbits,compress_buffer,&output_ptr); - write_instruction(coder,INSTR_DEFAULT,&output_ptr); - write_three_large(coder,has_large_ints,large_index,large_nbits,compress_buffer,&output_ptr); - has_large=0; - } - } - else - { + fprintf(stderr,"Sequence of one large and small integers (good compression).\n"); +#endif + /* Flush all large atoms but one! */ + if (has_large>1) + flush_large(coder,&has_large,has_large_ints,has_large-1,large_index,large_nbits,compress_buffer,&output_ptr); + write_instruction(coder,INSTR_DEFAULT,&output_ptr); + write_three_large(coder,has_large_ints,large_index,large_nbits,compress_buffer,&output_ptr); + has_large=0; + } + } + else + { #ifdef SHOWIT - fprintf(stderr,"Sequence of only small integers.\n"); -#endif - write_instruction(coder,INSTR_ONLY_SMALL,&output_ptr); - } - /* Base compress small integers using the current parameters. */ - nbits=magic_bits[small_index][runlength-1]; - /* The same base is used for the small changes. */ - small_idx[0]=small_index; - small_idx[1]=small_index; - small_idx[2]=small_index; - trajcoder_base_compress(encode_ints,runlength*3,small_idx,compress_buffer); + fprintf(stderr,"Sequence of only small integers.\n"); +#endif + write_instruction(coder,INSTR_ONLY_SMALL,&output_ptr); + } + /* Base compress small integers using the current parameters. */ + nbits=magic_bits[small_index][runlength-1]; + /* The same base is used for the small changes. */ + small_idx[0]=small_index; + small_idx[1]=small_index; + small_idx[2]=small_index; + trajcoder_base_compress(encode_ints,runlength*3,small_idx,compress_buffer); #ifdef SHOWIT - fprintf(stderr,"nbits=%d (%g)\n",nbits,nbits/(runlength*3.)); - nbits_sum+=nbits; - nvalues_sum+=runlength*3; - fprintf(stderr,"Runlength encoded small integers. runlength=%d\n",runlength); + fprintf(stderr,"nbits=%d (%g)\n",nbits,nbits/(runlength*3.)); + nbits_sum+=nbits; + nvalues_sum+=runlength*3; + fprintf(stderr,"Runlength encoded small integers. runlength=%d\n",runlength); #endif - /* write out base compressed small integers */ - Ptngc_writemanybits(coder,compress_buffer,nbits,&output_ptr); + /* write out base compressed small integers */ + Ptngc_writemanybits(coder,compress_buffer,nbits,&output_ptr); #ifdef SHOWIT - for (ienc=0; ienc=0) && (instr(unsigned int)normal_max) - normal_max=positive_int(normal[j]); - if (positive_int(swapped[j])>(unsigned int)swapped_max) - swapped_max=positive_int(swapped[j]); - } + { + if (positive_int(normal[j])>(unsigned int)normal_max) + normal_max=positive_int(normal[j]); + if (positive_int(swapped[j])>(unsigned int)swapped_max) + swapped_max=positive_int(swapped[j]); + } } if (normal_max==0) normal_max=1; @@ -202,8 +202,8 @@ static void allocate_enough_memory(unsigned int **ptr, int *nele, int *nele_allo } static void insert_value_in_array(unsigned int **ptr, int *nele, int *nele_alloc, - unsigned int value, - char *arrayname) + unsigned int value, + char *arrayname) { #ifndef SHOWIT (void)arrayname; @@ -237,21 +237,21 @@ static void swapdecide(struct xtc3_context *xtc3_context, int *input,int *swapat ((normalinstructions, - &xtc3_context->ninstr, - &xtc3_context->ninstr_alloc, - INSTR_FLIP,"instr"); + &xtc3_context->ninstr, + &xtc3_context->ninstr_alloc, + INSTR_FLIP,"instr"); } } @@ -277,11 +277,11 @@ static int is_quite_large(int *input, int small_index, int max_large_index) else { for (i=0; i<3; i++) - if (positive_int(input[i])>(unsigned int)Ptngc_magic(small_index+QUITE_LARGE)) - { - is=1; - break; - } + if (positive_int(input[i])>(unsigned int)Ptngc_magic(small_index+QUITE_LARGE)) + { + is=1; + break; + } } return is; } @@ -304,21 +304,21 @@ static void insert_batch(int *input_ptr, int ntriplets_left, int *prevcoord, int { int i; for (i=0; icurrent_large_type=xtc3_context->has_large_type[i]; if (xtc3_context->current_large_type==0) - instr=INSTR_LARGE_DIRECT; + instr=INSTR_LARGE_DIRECT; else if (xtc3_context->current_large_type==1) - instr=INSTR_LARGE_INTRA_DELTA; + instr=INSTR_LARGE_INTRA_DELTA; else - instr=INSTR_LARGE_INTER_DELTA; + instr=INSTR_LARGE_INTER_DELTA; insert_value_in_array(&xtc3_context->instructions, - &xtc3_context->ninstr, - &xtc3_context->ninstr_alloc, - instr,"instr"); + &xtc3_context->ninstr, + &xtc3_context->ninstr_alloc, + instr,"instr"); } } static void write_three_large(struct xtc3_context *xtc3_context, - int i) + int i) { int m; if (xtc3_context->current_large_type==0) { for (m=0; m<3; m++) - insert_value_in_array(&xtc3_context->large_direct, - &xtc3_context->nlargedir, - &xtc3_context->nlargedir_alloc, - xtc3_context->has_large_ints[i*3+m],"large direct"); + insert_value_in_array(&xtc3_context->large_direct, + &xtc3_context->nlargedir, + &xtc3_context->nlargedir_alloc, + xtc3_context->has_large_ints[i*3+m],"large direct"); } else if (xtc3_context->current_large_type==1) { for (m=0; m<3; m++) - insert_value_in_array(&xtc3_context->large_intra_delta, - &xtc3_context->nlargeintra, - &xtc3_context->nlargeintra_alloc, - xtc3_context->has_large_ints[i*3+m],"large intra"); + insert_value_in_array(&xtc3_context->large_intra_delta, + &xtc3_context->nlargeintra, + &xtc3_context->nlargeintra_alloc, + xtc3_context->has_large_ints[i*3+m],"large intra"); } else { for (m=0; m<3; m++) - insert_value_in_array(&xtc3_context->large_inter_delta, - &xtc3_context->nlargeinter, - &xtc3_context->nlargeinter_alloc, - xtc3_context->has_large_ints[i*3+m],"large inter"); + insert_value_in_array(&xtc3_context->large_inter_delta, + &xtc3_context->nlargeinter, + &xtc3_context->nlargeinter_alloc, + xtc3_context->has_large_ints[i*3+m],"large inter"); } } static void flush_large(struct xtc3_context *xtc3_context, - int n) /* How many to flush. */ + int n) /* How many to flush. */ { int i; i=0; @@ -409,48 +409,48 @@ static void flush_large(struct xtc3_context *xtc3_context, { int j,k; /* If the first large is of a different kind than the currently used we must - emit an "instruction" to change the large type. */ + emit an "instruction" to change the large type. */ large_instruction_change(xtc3_context,i); /* How many large of the same kind in a row? */ for (j=0; - (i+jhas_large_type[i+j]==xtc3_context->has_large_type[i]); - j++); + (i+jhas_large_type[i+j]==xtc3_context->has_large_type[i]); + j++); if (j<3) - { - for (k=0; kinstructions, - &xtc3_context->ninstr, - &xtc3_context->ninstr_alloc, - INSTR_ONLY_LARGE,"instr"); - write_three_large(xtc3_context,i+k); - } - } + { + for (k=0; kinstructions, + &xtc3_context->ninstr, + &xtc3_context->ninstr_alloc, + INSTR_ONLY_LARGE,"instr"); + write_three_large(xtc3_context,i+k); + } + } else - { - insert_value_in_array(&xtc3_context->instructions, - &xtc3_context->ninstr, - &xtc3_context->ninstr_alloc, - INSTR_LARGE_RLE,"instr"); - insert_value_in_array(&xtc3_context->rle, - &xtc3_context->nrle, - &xtc3_context->nrle_alloc, - (unsigned int)j,"rle (large)"); - for (k=0; kinstructions, + &xtc3_context->ninstr, + &xtc3_context->ninstr_alloc, + INSTR_LARGE_RLE,"instr"); + insert_value_in_array(&xtc3_context->rle, + &xtc3_context->nrle, + &xtc3_context->nrle_alloc, + (unsigned int)j,"rle (large)"); + for (k=0; khas_large-n)!=0) { int j; for (i=0; ihas_large-n; i++) - { - xtc3_context->has_large_type[i]=xtc3_context->has_large_type[i+n]; - for (j=0; j<3; j++) - xtc3_context->has_large_ints[i*3+j]=xtc3_context->has_large_ints[(i+n)*3+j]; - } + { + xtc3_context->has_large_type[i]=xtc3_context->has_large_type[i+n]; + for (j=0; j<3; j++) + xtc3_context->has_large_ints[i*3+j]=xtc3_context->has_large_ints[(i+n)*3+j]; + } } xtc3_context->has_large-=n; /* Number of remaining large atoms in buffer */ } @@ -467,7 +467,7 @@ static double compute_intlen(unsigned int *ints) } static void buffer_large(struct xtc3_context *xtc3_context, int *input, int inpdata, - int natoms, int intradelta_ok) + int natoms, int intradelta_ok) { unsigned int direct[3], intradelta[3]={0,}, interdelta[3]={0,}; double minlen; @@ -495,10 +495,10 @@ static void buffer_large(struct xtc3_context *xtc3_context, int *input, int inpd intradelta[2]=positive_int(input[inpdata+2]-input[inpdata-1]); thislen=compute_intlen(intradelta); if (thislen*TRESHOLD_INTRA_INTER_DIRECThas_large_type[xtc3_context->has_large]=best_type; @@ -581,18 +581,18 @@ static int base_bytes(unsigned int base, int n) for (i=0; i>(j*8))&0xFFU) - numbytes=i*4+j+1; + if ((largeint[i]>>(j*8))&0xFFU) + numbytes=i*4+j+1; return numbytes; } @@ -615,89 +615,89 @@ static void base_compress(unsigned int *data, int len, unsigned char *output, in int nvals=0; int basegiven=0; for (j=0; jbase) - base=data[k]; - basecheckvals++; - if (basecheckvals==MAXBASEVALS*BASEINTERVAL) - break; - } - /* The base is one larger than the largest values. */ - base++; - if (base<2) - base=2; - /* Store the base in the output. */ - output[nwrittenout++]=(unsigned char)(base&0xFFU); - output[nwrittenout++]=(unsigned char)((base>>8)&0xFFU); - output[nwrittenout++]=(unsigned char)((base>>16)&0xFFU); - output[nwrittenout++]=(unsigned char)((base>>24)&0xFFU); - basegiven=BASEINTERVAL; - /* How many bytes is needed to store MAXBASEVALS values using this base? */ - numbytes=base_bytes(base,MAXBASEVALS); - } - basegiven--; + if (basegiven==0) + { + base=0U; + /* Find the largest value for this particular coordinate. */ + for (k=i; kbase) + base=data[k]; + basecheckvals++; + if (basecheckvals==MAXBASEVALS*BASEINTERVAL) + break; + } + /* The base is one larger than the largest values. */ + base++; + if (base<2) + base=2; + /* Store the base in the output. */ + output[nwrittenout++]=(unsigned char)(base&0xFFU); + output[nwrittenout++]=(unsigned char)((base>>8)&0xFFU); + output[nwrittenout++]=(unsigned char)((base>>16)&0xFFU); + output[nwrittenout++]=(unsigned char)((base>>24)&0xFFU); + basegiven=BASEINTERVAL; + /* How many bytes is needed to store MAXBASEVALS values using this base? */ + numbytes=base_bytes(base,MAXBASEVALS); + } + basegiven--; #ifdef SHOWIT fprintf(stderr,"Base for %d is %u. I need %d bytes for %d values.\n",ixyz,base,numbytes,MAXBASEVALS); #endif } - if (nvals!=0) - { - Ptngc_largeint_mul(base,largeint,largeint_tmp,MAXBASEVALS+1); - for (j=0; j>(ibyte*8))&(0xFFU)); + for (j=0; j>(ibyte*8))&(0xFFU)); #ifdef SHOWIT - fprintf(stderr,"%02x",(unsigned int)output[nwrittenout-1]); + fprintf(stderr,"%02x",(unsigned int)output[nwrittenout-1]); #endif - } + } #ifdef SHOWIT - fprintf(stderr,"\n"); + fprintf(stderr,"\n"); #endif - nvals=0; - for (j=0; j>(ibyte*8))&(0xFFU)); - } - } + fprintf(stderr,"Base for %d is %u. I need %d bytes for %d values.\n",ixyz,base,numbytes,nvals); +#endif + for (j=0; j>(ibyte*8))&(0xFFU)); + } + } } *outlen=nwrittenout; } @@ -712,8 +712,8 @@ static void base_decompress(unsigned char *input, int len, unsigned int *output) if (maxbasevals>(int)MAXMAXBASEVALS) { fprintf(stderr,"Read a larger maxbasevals value from the file than I can handle. Fix" - " by increasing MAXMAXBASEVALS to at least %d. Although, this is" - " probably a bug in TRAJNG, since MAXMAXBASEVALS should already be insanely large enough.\n",maxbasevals); + " by increasing MAXMAXBASEVALS to at least %d. Although, this is" + " probably a bug in TRAJNG, since MAXMAXBASEVALS should already be insanely large enough.\n",maxbasevals); exit(EXIT_FAILURE); } input+=3; @@ -728,65 +728,65 @@ static void base_decompress(unsigned char *input, int len, unsigned int *output) fprintf(stderr,"Base for %d is %u. I need %d bytes for %d values.\n",ixyz,base,numbytes,maxbasevals); #endif while (nvals_left) - { - int n; - if (basegiven==0) - { - base=(unsigned int)(input[0])| - (((unsigned int)(input[1]))<<8)| - (((unsigned int)(input[2]))<<16)| - (((unsigned int)(input[3]))<<24); - input+=4; - basegiven=baseinterval; - /* How many bytes is needed to store maxbasevals values using this base? */ - numbytes=base_bytes(base,maxbasevals); - } - basegiven--; - if (nvals_leftnvals_left) - n=nvals_left; - for (i=n-1; i>=0; i--) - { - output[outvals+i*3]=Ptngc_largeint_div(base,largeint,largeint_tmp,maxbasevals+1); - for (j=0; jnvals_left) + n=nvals_left; + for (i=n-1; i>=0; i--) + { + output[outvals+i*3]=Ptngc_largeint_div(base,largeint,largeint_tmp,maxbasevals+1); + for (j=0; jxtc3_context.maxint[j]) - xtc3_context.maxint[j]=input[i*3+j]; - if (input[i*3+j]xtc3_context.maxint[j]) + xtc3_context.maxint[j]=input[i*3+j]; + if (input[i*3+j]intmax) - if (s0) - { - unsigned int delta[3], delta2[3]; - delta[0]=positive_int(input[inpdata+3]-input[inpdata-natoms*3+3]); - delta[1]=positive_int(input[inpdata+4]-input[inpdata-natoms*3+4]); - delta[2]=positive_int(input[inpdata+5]-input[inpdata-natoms*3+5]); - delta2[0]=positive_int(encode_ints[3]); - delta2[1]=positive_int(encode_ints[4]); - delta2[2]=positive_int(encode_ints[5]); + if (!no_swap) + { + /* If doing inter-frame coding results in smaller values we should not do any swapping either. */ + int frame=inpdata/(natoms*3); + if (frame>0) + { + unsigned int delta[3], delta2[3]; + delta[0]=positive_int(input[inpdata+3]-input[inpdata-natoms*3+3]); + delta[1]=positive_int(input[inpdata+4]-input[inpdata-natoms*3+4]); + delta[2]=positive_int(input[inpdata+5]-input[inpdata-natoms*3+5]); + delta2[0]=positive_int(encode_ints[3]); + delta2[1]=positive_int(encode_ints[4]); + delta2[2]=positive_int(encode_ints[5]); #ifdef SHOWIT - fprintf(stderr,"A1: inter delta: %u %u %u. intra delta=%u %u %u\n", - delta[0],delta[1],delta[2], - delta2[0],delta2[1],delta2[2]); -#endif - if (compute_intlen(delta)*TRESHOLD_INTER_INTRAlargest_required_base) - largest_required_base=encode_ints[ienc]; - /* Also compute what the largest base is for the current runlength setting! */ - largest_runlength_base=0; - for (ienc=0; (ienclargest_runlength_base) - largest_runlength_base=encode_ints[ienc]; - - largest_required_index=Ptngc_find_magic_index(largest_required_base); - largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); - - if (largest_required_indexntriplets_left) - new_runlength=ntriplets_left; - - /* We must at least try to get some small integers going. */ - if (new_runlength==0) - { - new_runlength=1; - new_small_index=small_index; - } - - iter_runlength=new_runlength; - iter_small_index=new_small_index; - - /* Iterate to find optimal encoding and runlength */ + fprintf(stderr,"Update batch due to large int.\n"); +#endif + if ((swapatoms) && (didswap)) + { + /* Keep swapped values. */ + for (i=0; i<2; i++) + for (ienc=0; ienc<3; ienc++) + encode_ints[i*3+ienc]=encode_ints[(i+1)*3+ienc]; + } + insert_batch(input+inpdata,ntriplets_left,prevcoord,encode_ints,min_runlength,&nencode); + } + /* Here we should only have differences for the atom coordinates. */ + /* Convert the ints to positive ints */ + for (ienc=0; ienclargest_required_base) + largest_required_base=encode_ints[ienc]; + /* Also compute what the largest base is for the current runlength setting! */ + largest_runlength_base=0; + for (ienc=0; (ienclargest_runlength_base) + largest_runlength_base=encode_ints[ienc]; + + largest_required_index=Ptngc_find_magic_index(largest_required_base); + largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); + + if (largest_required_indexntriplets_left) + new_runlength=ntriplets_left; + + /* We must at least try to get some small integers going. */ + if (new_runlength==0) + { + new_runlength=1; + new_small_index=small_index; + } + + iter_runlength=new_runlength; + iter_small_index=new_small_index; + + /* Iterate to find optimal encoding and runlength */ #ifdef SHOWIT - fprintf(stderr,"Entering iterative loop.\n"); - fflush(stderr); + fprintf(stderr,"Entering iterative loop.\n"); + fflush(stderr); #endif - do { - new_runlength=iter_runlength; - new_small_index=iter_small_index; + do { + new_runlength=iter_runlength; + new_small_index=iter_small_index; #ifdef SHOWIT - fprintf(stderr,"Test new_small_index=%d Base=%d\n",new_small_index,Ptngc_magic(new_small_index)); -#endif - /* What is the largest runlength - we can do with the currently - selected encoding? Also the max supported runlength is MAX_SMALL_RLE triplets! */ - for (ienc=0; iencnew_small_index) - break; - } - if (ienc/3>new_runlength) - { - iter_runlength=ienc/3; + fprintf(stderr,"Test new_small_index=%d Base=%d\n",new_small_index,Ptngc_magic(new_small_index)); +#endif + /* What is the largest runlength + we can do with the currently + selected encoding? Also the max supported runlength is MAX_SMALL_RLE triplets! */ + for (ienc=0; iencnew_small_index) + break; + } + if (ienc/3>new_runlength) + { + iter_runlength=ienc/3; #ifdef SHOWIT - fprintf(stderr,"I found a new possible runlength: %d\n",iter_runlength); -#endif - } - - /* How large encoding do we have to use? */ - largest_runlength_base=0; - for (ienc=0; ienclargest_runlength_base) - largest_runlength_base=encode_ints[ienc]; - largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); - if (largest_runlength_index!=new_small_index) - { - iter_small_index=largest_runlength_index; + fprintf(stderr,"I found a new possible runlength: %d\n",iter_runlength); +#endif + } + + /* How large encoding do we have to use? */ + largest_runlength_base=0; + for (ienc=0; ienclargest_runlength_base) + largest_runlength_base=encode_ints[ienc]; + largest_runlength_index=Ptngc_find_magic_index(largest_runlength_base); + if (largest_runlength_index!=new_small_index) + { + iter_small_index=largest_runlength_index; #ifdef SHOWIT - fprintf(stderr,"I found a new possible small index: %d Base=%d\n",iter_small_index,Ptngc_magic(iter_small_index)); + fprintf(stderr,"I found a new possible small index: %d Base=%d\n",iter_small_index,Ptngc_magic(iter_small_index)); #endif - } - } while ((new_runlength!=iter_runlength) || - (new_small_index!=iter_small_index)); + } + } while ((new_runlength!=iter_runlength) || + (new_small_index!=iter_small_index)); #ifdef SHOWIT - fprintf(stderr,"Exit iterative loop.\n"); - fflush(stderr); -#endif - - /* Verify that we got something good. We may have caught a - substantially larger atom. If so we should just bail - out and let the loop get on another lap. We may have a - minimum runlength though and then we have to fulfill - the request to write out these atoms! */ - rle_index_dep=0; - if (new_runlength<3) - rle_index_dep=IS_LARGE; - else if (new_runlength<6) - rle_index_dep=QUITE_LARGE; - if ((min_runlength) - || ((new_small_index0)) - { - for (i=0; i=2*new_runlength/3)) - { - /* Put all the values in large arrays, instead of the small array */ - if (new_runlength) - { - for (i=0; i0)) + { + for (i=0; i=2*new_runlength/3)) + { + /* Put all the values in large arrays, instead of the small array */ + if (new_runlength) + { + for (i=0; i=%g?\n",change,isum,(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)); + fprintf(stderr,"Tested decrease %d of index: %g>=%g?\n",change,isum,(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)); #endif - if (isum>(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)) - { + if (isum>(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)) + { #ifdef SHOWIT - fprintf(stderr,"Rejected decrease %d of index due to length of vector: %g>=%g\n",change,isum,(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)); -#endif - rejected=1; - change++; - } - } while ((change<0) && (rejected)); - if (change==0) - break; - } - } - - /* Always accept the new small indices here. */ - small_index=new_small_index; - /* If we have a new runlength emit it */ - if (runlength!=new_runlength) - { - runlength=new_runlength; - insert_value_in_array(&xtc3_context.instructions, - &xtc3_context.ninstr, - &xtc3_context.ninstr_alloc, - INSTR_SMALL_RUNLENGTH,"instr"); - insert_value_in_array(&xtc3_context.rle, - &xtc3_context.nrle, - &xtc3_context.nrle_alloc, - (unsigned int)runlength,"rle (small)"); - } + fprintf(stderr,"Rejected decrease %d of index due to length of vector: %g>=%g\n",change,isum,(double)Ptngc_magic(small_index+change)*(double)Ptngc_magic(small_index+change)); +#endif + rejected=1; + change++; + } + } while ((change<0) && (rejected)); + if (change==0) + break; + } + } + + /* Always accept the new small indices here. */ + small_index=new_small_index; + /* If we have a new runlength emit it */ + if (runlength!=new_runlength) + { + runlength=new_runlength; + insert_value_in_array(&xtc3_context.instructions, + &xtc3_context.ninstr, + &xtc3_context.ninstr_alloc, + INSTR_SMALL_RUNLENGTH,"instr"); + insert_value_in_array(&xtc3_context.rle, + &xtc3_context.nrle, + &xtc3_context.nrle_alloc, + (unsigned int)runlength,"rle (small)"); + } #ifdef SHOWIT - fprintf(stderr,"Current small index: %d Base=%d\n",small_index,Ptngc_magic(small_index)); -#endif - } - /* If we have a large previous integer we can combine it with a sequence of small ints. */ - if (xtc3_context.has_large) - { - /* If swapatoms is set to 1 but we did actually not - do any swapping, we must first write out the - large atom and then the small. If swapatoms is 1 - and we did swapping we can use the efficient - encoding. */ - if ((swapatoms) && (!didswap)) - { + fprintf(stderr,"Current small index: %d Base=%d\n",small_index,Ptngc_magic(small_index)); +#endif + } + /* If we have a large previous integer we can combine it with a sequence of small ints. */ + if (xtc3_context.has_large) + { + /* If swapatoms is set to 1 but we did actually not + do any swapping, we must first write out the + large atom and then the small. If swapatoms is 1 + and we did swapping we can use the efficient + encoding. */ + if ((swapatoms) && (!didswap)) + { #ifdef SHOWIT - fprintf(stderr,"Swapatoms was set to 1 but we did not do swapping!\n"); - fprintf(stderr,"Only one large integer.\n"); + fprintf(stderr,"Swapatoms was set to 1 but we did not do swapping!\n"); + fprintf(stderr,"Only one large integer.\n"); #endif - /* Flush all large atoms. */ - flush_large(&xtc3_context,xtc3_context.has_large); + /* Flush all large atoms. */ + flush_large(&xtc3_context,xtc3_context.has_large); #ifdef SHOWIT - fprintf(stderr,"Sequence of only small integers.\n"); + fprintf(stderr,"Sequence of only small integers.\n"); #endif - insert_value_in_array(&xtc3_context.instructions, - &xtc3_context.ninstr, - &xtc3_context.ninstr_alloc, - INSTR_ONLY_SMALL,"instr"); - } - else - { + insert_value_in_array(&xtc3_context.instructions, + &xtc3_context.ninstr, + &xtc3_context.ninstr_alloc, + INSTR_ONLY_SMALL,"instr"); + } + else + { #ifdef SHOWIT - fprintf(stderr,"Sequence of one large and small integers (good compression).\n"); -#endif - /* Flush all large atoms but one! */ - if (xtc3_context.has_large>1) - flush_large(&xtc3_context,xtc3_context.has_large-1); - - /* Here we must check if we should emit a large - type change instruction. */ - large_instruction_change(&xtc3_context,0); - - insert_value_in_array(&xtc3_context.instructions, - &xtc3_context.ninstr, - &xtc3_context.ninstr_alloc, - INSTR_DEFAULT,"instr"); - - write_three_large(&xtc3_context,0); - xtc3_context.has_large=0; - } - } - else - { + fprintf(stderr,"Sequence of one large and small integers (good compression).\n"); +#endif + /* Flush all large atoms but one! */ + if (xtc3_context.has_large>1) + flush_large(&xtc3_context,xtc3_context.has_large-1); + + /* Here we must check if we should emit a large + type change instruction. */ + large_instruction_change(&xtc3_context,0); + + insert_value_in_array(&xtc3_context.instructions, + &xtc3_context.ninstr, + &xtc3_context.ninstr_alloc, + INSTR_DEFAULT,"instr"); + + write_three_large(&xtc3_context,0); + xtc3_context.has_large=0; + } + } + else + { #ifdef SHOWIT - fprintf(stderr,"Sequence of only small integers.\n"); -#endif - insert_value_in_array(&xtc3_context.instructions, - &xtc3_context.ninstr, - &xtc3_context.ninstr_alloc, - INSTR_ONLY_SMALL,"instr"); - } - /* Insert the small integers into the small integer array. */ - for (ienc=0; ienc=5) - bwlzh_compress(xtc3_context.instructions,xtc3_context.ninstr,bwlzh_buf,&bwlzh_buf_len); + bwlzh_compress(xtc3_context.instructions,xtc3_context.ninstr,bwlzh_buf,&bwlzh_buf_len); else - bwlzh_compress_no_lz77(xtc3_context.instructions,xtc3_context.ninstr,bwlzh_buf,&bwlzh_buf_len); + bwlzh_compress_no_lz77(xtc3_context.instructions,xtc3_context.ninstr,bwlzh_buf,&bwlzh_buf_len); output_int(output,&outdata,(unsigned int)bwlzh_buf_len); memcpy(output+outdata,bwlzh_buf,bwlzh_buf_len); outdata+=bwlzh_buf_len; @@ -1474,9 +1474,9 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp { bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nrle)); if (speed>=5) - bwlzh_compress(xtc3_context.rle,xtc3_context.nrle,bwlzh_buf,&bwlzh_buf_len); + bwlzh_compress(xtc3_context.rle,xtc3_context.nrle,bwlzh_buf,&bwlzh_buf_len); else - bwlzh_compress_no_lz77(xtc3_context.rle,xtc3_context.nrle,bwlzh_buf,&bwlzh_buf_len); + bwlzh_compress_no_lz77(xtc3_context.rle,xtc3_context.nrle,bwlzh_buf,&bwlzh_buf_len); output_int(output,&outdata,(unsigned int)bwlzh_buf_len); memcpy(output+outdata,bwlzh_buf,bwlzh_buf_len); outdata+=bwlzh_buf_len; @@ -1491,18 +1491,18 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp if (xtc3_context.nlargedir) { if ((speed<=2) || ((speed<=5) && (!heuristic_bwlzh(xtc3_context.large_direct,xtc3_context.nlargedir)))) - { - bwlzh_buf=NULL; - bwlzh_buf_len=INT_MAX; - } + { + bwlzh_buf=NULL; + bwlzh_buf_len=INT_MAX; + } else - { - bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nlargedir)); - if (speed>=5) - bwlzh_compress(xtc3_context.large_direct,xtc3_context.nlargedir,bwlzh_buf,&bwlzh_buf_len); - else - bwlzh_compress_no_lz77(xtc3_context.large_direct,xtc3_context.nlargedir,bwlzh_buf,&bwlzh_buf_len); - } + { + bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nlargedir)); + if (speed>=5) + bwlzh_compress(xtc3_context.large_direct,xtc3_context.nlargedir,bwlzh_buf,&bwlzh_buf_len); + else + bwlzh_compress_no_lz77(xtc3_context.large_direct,xtc3_context.nlargedir,bwlzh_buf,&bwlzh_buf_len); + } /* If this can be written smaller using base compression we should do that. */ base_buf=warnmalloc((xtc3_context.nlargedir+3)*sizeof(int)); base_compress(xtc3_context.large_direct,xtc3_context.nlargedir,base_buf,&base_buf_len); @@ -1510,19 +1510,19 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp fprintf(stderr,"Large direct: Base len=%d. BWLZH len=%d\n",base_buf_len,bwlzh_buf_len); #endif if (base_buf_len=5) - bwlzh_compress(xtc3_context.large_intra_delta,xtc3_context.nlargeintra,bwlzh_buf,&bwlzh_buf_len); - else - bwlzh_compress_no_lz77(xtc3_context.large_intra_delta,xtc3_context.nlargeintra,bwlzh_buf,&bwlzh_buf_len); - } + { + bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nlargeintra)); + if (speed>=5) + bwlzh_compress(xtc3_context.large_intra_delta,xtc3_context.nlargeintra,bwlzh_buf,&bwlzh_buf_len); + else + bwlzh_compress_no_lz77(xtc3_context.large_intra_delta,xtc3_context.nlargeintra,bwlzh_buf,&bwlzh_buf_len); + } /* If this can be written smaller using base compression we should do that. */ base_buf=warnmalloc((xtc3_context.nlargeintra+3)*sizeof(int)); base_compress(xtc3_context.large_intra_delta,xtc3_context.nlargeintra,base_buf,&base_buf_len); @@ -1554,19 +1554,19 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp fprintf(stderr,"Large intra: Base len=%d. BWLZH len=%d\n",base_buf_len,bwlzh_buf_len); #endif if (base_buf_len=5) - bwlzh_compress(xtc3_context.large_inter_delta,xtc3_context.nlargeinter,bwlzh_buf,&bwlzh_buf_len); - else - bwlzh_compress_no_lz77(xtc3_context.large_inter_delta,xtc3_context.nlargeinter,bwlzh_buf,&bwlzh_buf_len); - } + { + bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nlargeinter)); + if (speed>=5) + bwlzh_compress(xtc3_context.large_inter_delta,xtc3_context.nlargeinter,bwlzh_buf,&bwlzh_buf_len); + else + bwlzh_compress_no_lz77(xtc3_context.large_inter_delta,xtc3_context.nlargeinter,bwlzh_buf,&bwlzh_buf_len); + } /* If this can be written smaller using base compression we should do that. */ base_buf=warnmalloc((xtc3_context.nlargeinter+3)*sizeof(int)); base_compress(xtc3_context.large_inter_delta,xtc3_context.nlargeinter,base_buf,&base_buf_len); @@ -1598,19 +1598,19 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp fprintf(stderr,"Large inter: Base len=%d. BWLZH len=%d\n",base_buf_len,bwlzh_buf_len); #endif if (base_buf_len=5) - bwlzh_compress(xtc3_context.smallintra,xtc3_context.nsmallintra,bwlzh_buf,&bwlzh_buf_len); - else - bwlzh_compress_no_lz77(xtc3_context.smallintra,xtc3_context.nsmallintra,bwlzh_buf,&bwlzh_buf_len); - } + { + bwlzh_buf=warnmalloc(bwlzh_get_buflen(xtc3_context.nsmallintra)); + if (speed>=5) + bwlzh_compress(xtc3_context.smallintra,xtc3_context.nsmallintra,bwlzh_buf,&bwlzh_buf_len); + else + bwlzh_compress_no_lz77(xtc3_context.smallintra,xtc3_context.nsmallintra,bwlzh_buf,&bwlzh_buf_len); + } /* If this can be written smaller using base compression we should do that. */ base_buf=warnmalloc((xtc3_context.nsmallintra+3)*sizeof(int)); base_compress(xtc3_context.smallintra,xtc3_context.nsmallintra,base_buf,&base_buf_len); @@ -1642,19 +1642,19 @@ unsigned char *Ptngc_pack_array_xtc3(int *input, int *length, int natoms, int sp fprintf(stderr,"Small intra: Base len=%d. BWLZH len=%d\n",base_buf_len,bwlzh_buf_len); #endif if (base_buf_lenlarge_direct) @@ -1717,11 +1717,11 @@ static void unpack_one_large(struct xtc3_context *xtc3_context, else if (xtc3_context->large_inter_delta) { large_ints[0]=unpositive_int(xtc3_context->large_inter_delta[(*ilargeinter)]) - +output[outdata-natoms*3+didswap*3]; + +output[outdata-natoms*3+didswap*3]; large_ints[1]=unpositive_int(xtc3_context->large_inter_delta[(*ilargeinter)+1]) - +output[outdata-natoms*3+1+didswap*3]; + +output[outdata-natoms*3+1+didswap*3]; large_ints[2]=unpositive_int(xtc3_context->large_inter_delta[(*ilargeinter)+2]) - +output[outdata-natoms*3+2+didswap*3]; + +output[outdata-natoms*3+2+didswap*3]; (*ilargeinter)+=3; } prevcoord[0]=large_ints[0]; @@ -1760,78 +1760,78 @@ int Ptngc_unpack_array_xtc3(unsigned char *packed,int *output, int length, int n for (i=0; i<3; i++) { minint[i]=unpositive_int((int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24))); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24))); ptr+=4; } xtc3_context.ninstr=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.ninstr) decompress_bwlzh_block(&ptr,xtc3_context.ninstr,&xtc3_context.instructions); xtc3_context.nrle=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.nrle) decompress_bwlzh_block(&ptr,xtc3_context.nrle,&xtc3_context.rle); xtc3_context.nlargedir=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.nlargedir) { if (*ptr++==1) - decompress_bwlzh_block(&ptr,xtc3_context.nlargedir,&xtc3_context.large_direct); + decompress_bwlzh_block(&ptr,xtc3_context.nlargedir,&xtc3_context.large_direct); else - decompress_base_block(&ptr,xtc3_context.nlargedir,&xtc3_context.large_direct); + decompress_base_block(&ptr,xtc3_context.nlargedir,&xtc3_context.large_direct); } xtc3_context.nlargeintra=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.nlargeintra) { if (*ptr++==1) - decompress_bwlzh_block(&ptr,xtc3_context.nlargeintra,&xtc3_context.large_intra_delta); + decompress_bwlzh_block(&ptr,xtc3_context.nlargeintra,&xtc3_context.large_intra_delta); else - decompress_base_block(&ptr,xtc3_context.nlargeintra,&xtc3_context.large_intra_delta); + decompress_base_block(&ptr,xtc3_context.nlargeintra,&xtc3_context.large_intra_delta); } xtc3_context.nlargeinter=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.nlargeinter) { if (*ptr++==1) - decompress_bwlzh_block(&ptr,xtc3_context.nlargeinter,&xtc3_context.large_inter_delta); + decompress_bwlzh_block(&ptr,xtc3_context.nlargeinter,&xtc3_context.large_inter_delta); else - decompress_base_block(&ptr,xtc3_context.nlargeinter,&xtc3_context.large_inter_delta); + decompress_base_block(&ptr,xtc3_context.nlargeinter,&xtc3_context.large_inter_delta); } xtc3_context.nsmallintra=(int)(((unsigned int)ptr[0]) | - (((unsigned int)ptr[1])<<8) | - (((unsigned int)ptr[2])<<16) | - (((unsigned int)ptr[3])<<24)); + (((unsigned int)ptr[1])<<8) | + (((unsigned int)ptr[2])<<16) | + (((unsigned int)ptr[3])<<24)); ptr+=4; if (xtc3_context.nsmallintra) { if (*ptr++==1) - decompress_bwlzh_block(&ptr,xtc3_context.nsmallintra,&xtc3_context.smallintra); + decompress_bwlzh_block(&ptr,xtc3_context.nsmallintra,&xtc3_context.smallintra); else - decompress_base_block(&ptr,xtc3_context.nsmallintra,&xtc3_context.smallintra); + decompress_base_block(&ptr,xtc3_context.nsmallintra,&xtc3_context.smallintra); } /* Initial prevcoord is the minimum integers. */ @@ -1849,103 +1849,103 @@ int Ptngc_unpack_array_xtc3(unsigned char *packed,int *output, int length, int n fprintf(stderr,"ntriplets left=%d\n",ntriplets_left); #endif if ((instr==INSTR_DEFAULT) /* large+small */ - || (instr==INSTR_ONLY_LARGE) /* only large */ - || (instr==INSTR_ONLY_SMALL)) /* only small */ - { - if (instr!=INSTR_ONLY_SMALL) - { - int didswap=0; - if ((instr==INSTR_DEFAULT) && (swapatoms)) - didswap=1; - unpack_one_large(&xtc3_context,&ilargedir, &ilargeintra, &ilargeinter, - prevcoord, minint, output, outdata, didswap, - natoms, current_large_type); - ntriplets_left--; - outdata+=3; - } - if (instr!=INSTR_ONLY_LARGE) - { - for (i=0; iinput_file) - + file_pos = (int64_t)ftell(tng_data->input_file) - (block->block_contents_size + block->header_contents_size); if(hash_mode == TNG_USE_HASH) @@ -3863,7 +3863,7 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, { dest = tng_compress_pos_float_find_algo(start_pos, (int)n_particles, (int)algo_find_n_frames, - 0.001, 0, + (float)0.001, 0, tng_data-> compress_algo_pos, &new_len); @@ -3871,7 +3871,7 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, if(algo_find_n_frames < n_frames) { dest = tng_compress_pos_float(start_pos, (int)n_particles, - (int)n_frames, 0.001, 0, + (int)n_frames, (float)0.001, 0, tng_data->compress_algo_pos, &new_len); } @@ -3898,7 +3898,7 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, if(type == TNG_FLOAT_DATA) { dest = tng_compress_pos_float(start_pos, (int)n_particles, - (int)n_frames, 0.001, 0, + (int)n_frames, (float)0.001, 0, tng_data->compress_algo_pos, &new_len); } else @@ -3930,14 +3930,14 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, { dest = tng_compress_vel_float_find_algo(start_pos, (int)n_particles, (int)algo_find_n_frames, - 0.001, 0, + (float)0.001, 0, tng_data-> compress_algo_vel, &new_len); if(algo_find_n_frames < n_frames) { dest = tng_compress_vel_float(start_pos, (int)n_particles, - (int)n_frames, 0.001, 0, + (int)n_frames, (float)0.001, 0, tng_data->compress_algo_pos, &new_len); } @@ -3964,7 +3964,7 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, if(type == TNG_FLOAT_DATA) { dest = tng_compress_vel_float(start_pos, (int)n_particles, - (int)n_frames, 0.001, 0, + (int)n_frames, (float)0.001, 0, tng_data-> compress_algo_vel, &new_len); @@ -3985,7 +3985,7 @@ static tng_function_status tng_compress(tng_trajectory_t tng_data, return(TNG_FAILURE); } - offset = (char *)start_pos - block->block_contents; + offset = (unsigned long)((char *)start_pos - block->block_contents); block->block_contents_size = new_len + offset; @@ -4022,7 +4022,8 @@ static tng_function_status tng_uncompress(tng_trajectory_t tng_data, char *temp; double *d_dest = 0; float *f_dest = 0; - int offset, result; + unsigned long offset; + int result; (void)tng_data; if(block->id != TNG_TRAJ_POSITIONS && @@ -4068,10 +4069,10 @@ static tng_function_status tng_uncompress(tng_trajectory_t tng_data, return(TNG_FAILURE); } - offset = (char *)start_pos - (char *)block->block_contents; + offset = (unsigned long)((char *)start_pos - (char *)block->block_contents); - block->block_contents_size = uncompressed_len + offset; + block->block_contents_size = (int64_t)(uncompressed_len + offset); temp = realloc(block->block_contents, uncompressed_len + offset); if(!temp) @@ -4120,7 +4121,7 @@ static tng_function_status tng_gzip_compress(tng_trajectory_t tng_data, { Bytef *dest; char *temp; - uLong max_len, stat, offset; + unsigned long max_len, stat, offset; (void)tng_data; max_len = compressBound(len); @@ -4593,7 +4594,7 @@ static tng_function_status tng_particle_data_read if(codec_id != TNG_UNCOMPRESSED) { - data_size = n_frames_div * size * n_particles * n_values; + data_size = (unsigned long)(n_frames_div * size * n_particles * n_values); switch(codec_id) { case TNG_XTC_COMPRESSION: @@ -4739,10 +4740,10 @@ static tng_function_status tng_particle_data_block_write const char hash_mode) { int64_t n_particles, num_first_particle, n_frames, stride_length; - int64_t frame_step; + int64_t frame_step, data_start_pos; int64_t i, j, k; - int offset = 0, size, data_start_pos; - unsigned int len; + int offset = 0, size; + size_t len; char dependency, temp, *temp_name; double multiplier; char ***first_dim_values, **second_dim_values; @@ -4803,7 +4804,7 @@ static tng_function_status tng_particle_data_block_write temp_name = realloc(block->name, len); if(!temp_name) { - printf("Cannot allocate memory (%d bytes). %s: %d\n", len, + printf("Cannot allocate memory (%ld bytes). %s: %d\n", len, __FILE__, __LINE__); free(block->name); block->name = 0; @@ -5096,7 +5097,7 @@ static tng_function_status tng_particle_data_block_write { for(i = offset; i < block->block_contents_size; i+=size) { - *(float *)(block->block_contents + i) *= multiplier; + *(float *)(block->block_contents + i) *= (float)multiplier; if(tng_data->input_endianness_swap_func_32 && tng_data->input_endianness_swap_func_32(tng_data, (int32_t *)(block->block_contents + i)) @@ -5183,7 +5184,7 @@ static tng_function_status tng_particle_data_block_write data->codec_id = TNG_UNCOMPRESSED; break; case TNG_TNG_COMPRESSION: - printf("TEST: frame_step: %"PRId64", n_particles: %"PRId64", block_contents_size: %"PRId64", start_pos: %d\n", + printf("TEST: frame_step: %"PRId64", n_particles: %"PRId64", block_contents_size: %"PRId64", start_pos: %"PRId64"\n", frame_step, n_particles, block->block_contents_size, data_start_pos); stat = tng_compress(tng_data, block, frame_step, n_particles, data->datatype, @@ -7059,8 +7060,7 @@ tng_function_status DECLSPECDLLEXPORT tng_molecule_cnt_get tng_molecule_t molecule, int64_t *cnt) { - int index = -1; - int64_t i; + int64_t i, index = -1; for(i = tng_data->n_molecules; i--;) { @@ -7084,8 +7084,7 @@ tng_function_status DECLSPECDLLEXPORT tng_molecule_cnt_set tng_molecule_t molecule, const int64_t cnt) { - int index = -1; - int64_t i, old_cnt; + int64_t i, old_cnt, index = -1; for(i = tng_data->n_molecules; i--;) { @@ -7307,8 +7306,7 @@ tng_function_status DECLSPECDLLEXPORT tng_chain_residue_w_id_add const int64_t id, tng_residue_t *residue) { - int64_t i; - int curr_index; + int64_t i, curr_index; tng_residue_t new_residues, temp_residue, last_residue; tng_molecule_t molecule = chain->molecule; tng_function_status stat = TNG_SUCCESS; @@ -9585,7 +9583,7 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_set_nr_find fseek(tng_data->input_file, (long)file_pos, SEEK_SET); - tng_data->current_trajectory_frame_set_input_file_pos = file_pos; + tng_data->current_trajectory_frame_set_input_file_pos = (long)file_pos; /* Read block headers first to see what block is found. */ stat = tng_block_header_read(tng_data, block); if(stat == TNG_CRITICAL || block->id != TNG_TRAJECTORY_FRAME_SET) @@ -9912,7 +9910,7 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_set_of_frame_find fseek(tng_data->input_file, (long)file_pos, SEEK_SET); - tng_data->current_trajectory_frame_set_input_file_pos = file_pos; + tng_data->current_trajectory_frame_set_input_file_pos = (long)file_pos; /* Read block headers first to see what block is found. */ stat = tng_block_header_read(tng_data, block); if(stat == TNG_CRITICAL || block->id != TNG_TRAJECTORY_FRAME_SET) @@ -10435,7 +10433,7 @@ tng_function_status DECLSPECDLLEXPORT tng_block_read_next(tng_trajectory_t tng_d tng_function_status DECLSPECDLLEXPORT tng_frame_set_read_next(tng_trajectory_t tng_data, const char hash_mode) { - long int file_pos; + long file_pos; tng_gen_block_t block; tng_function_status stat = TNG_SUCCESS; @@ -10444,11 +10442,11 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_set_read_next(tng_trajectory_t t return(TNG_CRITICAL); } - file_pos = (long int)tng_data->current_trajectory_frame_set.next_frame_set_file_pos; + file_pos = (long)tng_data->current_trajectory_frame_set.next_frame_set_file_pos; if(file_pos < 0 && tng_data->current_trajectory_frame_set_input_file_pos <= 0) { - file_pos = tng_data->first_trajectory_frame_set_input_file_pos; + file_pos = (long)tng_data->first_trajectory_frame_set_input_file_pos; } if(file_pos > 0) @@ -10535,8 +10533,9 @@ tng_function_status tng_frame_set_write(tng_trajectory_t tng_data, tng_function_status stat; tng_data->current_trajectory_frame_set_output_file_pos = - tng_data->last_trajectory_frame_set_output_file_pos = ftell(tng_data->output_file); + tng_data->last_trajectory_frame_set_output_file_pos = + tng_data->current_trajectory_frame_set_output_file_pos; if(tng_data->current_trajectory_frame_set_output_file_pos <= 0) { @@ -10613,7 +10612,7 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_set_new const int64_t first_frame, const int64_t n_frames) { - int i; + int64_t i; tng_gen_block_t block; tng_trajectory_frame_set_t frame_set; tng_particle_mapping_t mapping; @@ -10987,7 +10986,8 @@ tng_function_status DECLSPECDLLEXPORT tng_particle_data_block_add const int64_t codec_id, void *new_data) { - int i, j, k, size, len; + int i, size, len; + int64_t j, k; int64_t tot_n_particles, n_frames_div; char ***first_dim_values, **second_dim_values; tng_trajectory_frame_set_t frame_set; @@ -11140,8 +11140,9 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_data_write { int64_t header_pos, file_pos; int64_t output_file_len, n_values_per_frame, size, contents_size; - int64_t header_size, temp_first, temp_last, temp_current; + int64_t header_size, temp_first, temp_last; int64_t i, last_frame; + long temp_current; tng_gen_block_t block; tng_trajectory_frame_set_t frame_set; FILE *temp = tng_data->input_file; @@ -11553,9 +11554,10 @@ tng_function_status DECLSPECDLLEXPORT tng_frame_particle_data_write { int64_t header_pos, file_pos, tot_n_particles; int64_t output_file_len, n_values_per_frame, size, contents_size; - int64_t header_size, temp_first, temp_last, temp_current; + int64_t header_size, temp_first, temp_last; int64_t mapping_block_end_pos, num_first_particle, block_n_particles; int64_t i, last_frame; + long temp_current; tng_gen_block_t block; tng_trajectory_frame_set_t frame_set; FILE *temp = tng_data->input_file; @@ -12134,7 +12136,7 @@ tng_function_status DECLSPECDLLEXPORT tng_data_values_free const int64_t n_values_per_frame, const char type) { - int i, j; + int64_t i, j; (void)tng_data; if(values) @@ -12244,7 +12246,7 @@ tng_function_status DECLSPECDLLEXPORT tng_particle_data_values_free const int64_t n_values_per_frame, const char type) { - int i, j, k; + int64_t i, j, k; (void)tng_data; if(values) @@ -12289,8 +12291,8 @@ tng_function_status DECLSPECDLLEXPORT tng_data_get int64_t *n_values_per_frame, char *type) { - int64_t file_pos, block_index; - int i, j, len, size; + int64_t i, j, file_pos, block_index; + int len, size; tng_non_particle_data_t data; tng_trajectory_frame_set_t frame_set = &tng_data->current_trajectory_frame_set; @@ -12531,7 +12533,8 @@ tng_function_status DECLSPECDLLEXPORT tng_data_interval_get { int64_t i, j, n_frames, file_pos, current_frame_pos, first_frame; int64_t block_index; - int len, size; + int size; + size_t len; tng_non_particle_data_t data; tng_trajectory_frame_set_t frame_set; tng_gen_block_t block; @@ -13267,7 +13270,8 @@ tng_function_status DECLSPECDLLEXPORT tng_particle_data_interval_get { int64_t i, j, k, mapping, n_frames, file_pos, current_frame_pos, i_step; int64_t first_frame, block_index; - int len, size; + int size; + size_t len; tng_particle_data_t data; tng_trajectory_frame_set_t frame_set; tng_gen_block_t block; -- cgit v0.10.1