#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
#include <limits.h>
#include <math.h>

#define MAX(a,b)            (((a)>(b)) ? (a) : (b))

/* Initialize source buffer of integer for integer->integer and conversion test.
 * This algorithm is mainly to avoid any casting and comparison between source and destination types
 * for compiler, because we're testing conversions. */
#define INIT_INTEGER(TYPE, SRC_MAX, SRC_MIN, SRC_SIZE, DST_SIZE, SRC_PREC, BUF, NELMTS)  \
{                                                                                               \
    unsigned char *buf_p;                                                             \
    unsigned int n;                                                                             \
    TYPE value1 = 1;                                                                            \
    TYPE value2 = 0;                                                                            \
                                                                                                \
    /* Allocate buffers */                                                                      \
    NELMTS=(SRC_PREC-1)*3+1;                                                                    \
    BUF = (unsigned char*)malloc(NELMTS*MAX(SRC_SIZE, DST_SIZE));                       \
    memset(BUF, 0, NELMTS*MAX(SRC_SIZE, DST_SIZE));                                           \
                                                                                                \
    buf_p = BUF;                                                                                \
                                                                                                \
    /*positive values, ascending order. VALUE1 starts from 00000001, to 00000010, until 10000000*/ \
    /*VALUE2 ascends from 00000000, to 00000011, 00000111,...,  until 11111111.*/               \
    for(n=0; n<SRC_PREC-1; n++) {                                                               \
        if(value1<=SRC_MAX && value1>=SRC_MIN) {                                                \
            memcpy(buf_p, &value1, SRC_SIZE);                                                   \
            buf_p += SRC_SIZE;                                                                  \
        }                                                                                       \
        if(value2<=SRC_MAX && value2>=SRC_MIN) {                                                \
            memcpy(buf_p, &value2, SRC_SIZE);                                                   \
            buf_p += SRC_SIZE;                                                                  \
        }                                                                                       \
                                                                                                \
        value1 <<= 1;                                                                           \
        value2 = (value1 - 1) | value1;                                                         \
    }                                                                                           \
                                                                                                \
    /* negative values for signed; descending positive values for unsigned */                   \
    /* VALUE2 descends from 11111111 to 11111110, 11111100, ..., until 10000000. */             \
    for(n=0; n<SRC_PREC; n++) {                                                                 \
        if(value2<=SRC_MAX && value2>=SRC_MIN) {                                                \
            memcpy(buf_p, &value2, SRC_SIZE);                                                   \
            buf_p += SRC_SIZE;                                                                  \
        }                                                                                       \
        value2 <<= 1;                                                                           \
    }                                                                                           \
}

/* Convert from short to integer using the compiler. */
conv_short_int( size_t nelmts, void *buf)
{
      size_t elmtno; 
      uint8_t *src, *s; 
      uint8_t *dst, *d; 
      long long aligned; 
      ssize_t s_stride = 2;
      ssize_t d_stride = 4; 

      /* In-place conversion, from the last to the first element. */
      if(d_stride>s_stride) { 
	  src = (uint8_t*)buf+(nelmts-1)*s_stride; 
	  dst = (uint8_t*)buf+(nelmts-1)*d_stride; 
	  s_stride = -s_stride; 
	  d_stride = -d_stride; 
      } 

      for (elmtno=0; elmtno<nelmts; elmtno++) { 
	  s = (uint8_t*)&aligned; 
	  memcpy((char*)(&aligned),(const char*)(src),sizeof(short)); 
	  d = (uint8_t*)&aligned; 
 
          /* This line causes the trouble. */
	  *((int*)d) = (int)(*((short*)s));

	  memcpy((char*)(dst),(const char*)(&aligned),sizeof(int));
	  src += s_stride; 
	  dst += d_stride;
      } 

}

/* This C program simulates what the data conversion of the HDF5 software does.
 * It converts some data from short to int in place.  When I use -O2 or -O3 to 
 * compile it on linux, I got incorrect result.  -O or -O1 works.  I used 
 * Purify to debug it.  It pointed to line 78 as an error but didn't provide 
 * much information.
 */
int main(void)
{
    size_t	        nelmts=0;		/*num values per test	*/
    size_t		src_size, dst_size;	/*type sizes		*/
    unsigned char	*buf=NULL;		/*buffer for conversion	*/
    unsigned char	*p=NULL;		/*buffer for conversion	*/
    size_t		i;			/*counters		*/
    unsigned char	src_bits[32];		/*src value in LE order	*/
    unsigned char	dst_bits[32];		/*dest value in LE order*/
    size_t		src_nbits;		/*source length in bits	*/
    size_t		dst_nbits;		/*dst length in bits	*/

    src_size = sizeof(short);
    dst_size = sizeof(int);
    src_nbits = 8*sizeof(short);

    /* Allocate the buffer big enough for both source and destination data and 
     * initialize it with source data. */    
    INIT_INTEGER(short, SHRT_MAX, SHRT_MIN, src_size, dst_size, src_nbits, buf, nelmts);

    /* Print out the data before the conversion */
    p = buf;
    printf("before conversion:\n");
    for(i=0; i<nelmts; i++) {
       printf("%d. %d\n", i, *((short*)p)); 
       p += src_size;
    }

    /* The actual conversion */
    conv_short_int(nelmts, buf);

    /* Print out the data after the conversion */
    p = buf;
    printf("after conversion:\n");
    for(i=0; i<nelmts; i++) {
       printf("%d. %d\n", i, *((int*)p)); 
       p += dst_size;
    }

    if (buf) free(buf);

    return 0;
}