/*******************************************************************/
/*                                                                 */
/*  sd_create.c                                                    */
/*        This program creates an HDF file, file1.hdf, which          */
/*        contains an SDS of size 10 x 5.                          */
/*                                                                 */
/*******************************************************************/

#include "mfhdf.h"

#define DIM0 4
#define DIM1 3
#define DIM2 2
#define XX 5
#define YY 5

#define CHECKFAIL(status, mesg) {\
	if (status == -1) \
	    fprintf(stderr," %s\n", mesg);}

main( ) 
{
	int32 sd_id, sds1_id, sds2_id, sds_id3, sds_idx, sds_idx1, istat;
	int32 dims[3], start[3], start_w[3], rank=3;
	int32 edges[3], edges1[3], edges2[1];
	int16 data1_slab3_4[DIM0][DIM1][DIM2], 
		data1_slab1[DIM0][DIM1][DIM2], 
		data1_slab2[DIM0][DIM1][DIM2], 
		data1_slab234[DIM0][DIM1][DIM2], 
		outdata[DIM0][DIM1][DIM2],
		data2_slab1[DIM0][DIM1],
		data2_slab2[DIM0][DIM1],
		outdata2[DIM0][DIM1];
	char name1[10], name2[10];
	int32 dtype, nattrs, dimsizes[3];
	intn ii, jj, kk;
	int16 fillval = 99;

	/* Create and open the file and initiate the SD interface. */
	sd_id = SDstart("hdffile_unlim.hdf", DFACC_CREATE);
	if (sd_id == FAIL) {
		fprintf(stderr,"SDstart failed.\n");
		HEprint(stdout,0);
		exit (-1);
	}
	else
		fprintf(stderr,"\n... hdffile_unlim.hdf created\n");

	/* Define the rank and dimensions of the data set to be created. */
	rank = 3;
	dims[0] = SD_UNLIMITED;
	dims[1] = DIM1;
	dims[2] = DIM2;

	/* Create the array data sets. */
	sds1_id = SDcreate(sd_id, "Dataset 1", DFNT_INT16, 3, dims);
	istat = SDsetfillvalue(sds1_id, (VOIDP)&fillval);
	CHECKFAIL(istat, " SDsetfillvalue Dataset 1 failed \n");

	sds2_id = SDcreate(sd_id, "Dataset 2", DFNT_INT16, 2, dims);
	istat = SDsetfillvalue(sds2_id, (VOIDP)&fillval);
	CHECKFAIL(istat, " SDsetfillvalue Dataset 2 failed \n");

#if 0
	istat = SDendaccess(sds_id);

        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds_id = SDselect (sd_id, sds_idx);
#endif

	/* Fill the data arrays with values. */
	for (ii = 0; ii < 1; ii++)
	    for (jj = 0; jj < DIM1; jj++)
		for (kk = 0; kk < DIM2; kk++) {
		    data1_slab1[ii][jj][kk] =  kk + jj + 100;
		    data1_slab2[ii][jj][kk] =  kk + jj + (ii+2)*100;
		    }

	for (ii = 0; ii < 1; ii++)
	    for (jj = 0; jj < DIM1; jj++)
		    data2_slab1[ii][jj] = jj + 200;

/***** write first 1x3x2 slab *****/

	/* Define the location, pattern, and size of the data set */
	for (ii = 0; ii < rank; ii++) {
		start[ii] = 0;
		start_w[ii] = 0;
	}

	edges1[0] = 1;
	edges1[1] = DIM1;
	edges1[2] = DIM2;
        
	/* write the first time to Dataset 1 1x3x2*/
	istat = SDwritedata(sds1_id, start_w, NULL, edges1, (VOIDP)data1_slab1); 
	CHECKFAIL(istat, " SDwritedata Dataset 1 first time failed \n");

	/* write the first time to Dataset 2 1x3 */
	istat = SDwritedata(sds2_id, start_w, NULL, edges1, (VOIDP)data2_slab1); 
	CHECKFAIL(istat, " SDwritedata Dataset 2 first time failed \n");
	istat = SDendaccess(sds1_id);
	istat = SDendaccess(sds2_id);

	fprintf(stderr, "\n**** After writing first slab ****\n");

	/* Read and display Dataset 1 */
        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds1_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
	fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, 
		rank, dimsizes[0], dimsizes[1], dimsizes[2]);
	HDmemset(&outdata, 0, sizeof(outdata));
	for (ii = 0; ii < rank; ii++)
	    edges[ii] = dimsizes[ii];
        istat = SDreaddata(sds1_id, start, NULL, edges, (VOIDP) outdata);
	if (istat == -1) fprintf(stderr," SDreaddata failed \n");
        fprintf(stderr, "\nData read from Dataset 1 after first write: \n");
	for (ii = 0; ii < dimsizes[0]; ii++) {
	    for (jj = 0; jj < DIM1; jj++) {
		for (kk = 0; kk < DIM2; kk++) {
 		fprintf(stderr, " %d ", outdata[ii][jj][kk]);
		}
		fprintf(stderr, "\n");
 		}
 	    fprintf(stderr, "\n");
	}

	/* Read and display Dataset 2 */
        sds_idx = SDnametoindex (sd_id, "Dataset 2");
        sds2_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds2_id, name2, &rank, dimsizes, &dtype, &nattrs);
	fprintf(stderr, "dataset %s has %d dimensions: %d,%d\n", name2, 
		rank, dimsizes[0], dimsizes[1]);
	HDmemset(&outdata2, 0, sizeof(outdata2));
        istat = SDreaddata(sds2_id, start, NULL, edges, (VOIDP) outdata2);
	CHECKFAIL(istat, " SDreaddata failed \n");
        fprintf(stderr, "\nData read from Dataset 2 after first write: \n");
	for (ii = 0; ii < dimsizes[0]; ii++) {
	    for (jj = 0; jj < DIM1; jj++) {
 		fprintf(stderr, " %d ", outdata2[ii][jj]);
		}
 	    fprintf(stderr, "\n");
	}

	/* Close datasets */
	istat = SDendaccess(sds1_id);
	istat = SDendaccess(sds2_id);

	/* Display sizes before closing file */
        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds1_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, rank, dimsizes[0], dimsizes[1], dimsizes[2]);
	istat = SDendaccess(sds1_id);

        sds_idx = SDnametoindex (sd_id, "Dataset 2");
        sds2_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds2_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d\n", name1, rank, dimsizes[0], dimsizes[1]);
	istat = SDendaccess(sds2_id);

	/* close the file */
	istat = SDend(sd_id);

	/* Reopen the file for writing. */
	sd_id = SDstart("hdffile_unlim.hdf", DFACC_RDWR);
	if (sd_id == FAIL) {
		fprintf(stderr,"SDstart failed.\n");
		HEprint(stdout,0);
		exit (-1);
	}
	else
		fprintf(stderr,"\n... hdffile_unlim.hdf opened first time\n");

/***** write second 1x3x2 slab *****/

	/* Display dimension sizes after re-opening file */
        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds1_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, rank, dimsizes[0], dimsizes[1], dimsizes[2]);

	/* Define the location, pattern, and size of the data set */
	start_w[0] = 1;

	edges1[0] = 3;
	edges1[1] = DIM1;
	edges1[2] = DIM2;
        
	/* Fill the data arrays with values. */
	for (ii = 0; ii < 3; ii++)
	    for (jj = 0; jj < DIM1; jj++)
		for (kk = 0; kk < DIM2; kk++) {
		    data1_slab234[ii][jj][kk] =  kk + (jj+2)*10 + (ii+2)*100;
		    }

	/* write the second time to Dataset 1 */
	istat = SDwritedata(sds1_id, start_w, NULL, edges1, (VOIDP)data1_slab234); 
	CHECKFAIL(istat, " SDwritedata Dataset 1 first time failed \n");
 fprintf(stderr, "\n**** After writing second slab ****\n");

	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, rank, dimsizes[0], dimsizes[1], dimsizes[2]);

	HDmemset(&outdata, 0, sizeof(outdata));
	edges1[0] = 4; /* already wrote 4 slabs */
        istat = SDreaddata(sds1_id, start, NULL, edges1, (VOIDP) outdata);
	CHECKFAIL(istat, " SDreaddata failed \n");
        fprintf(stderr, "\nData read from Dataset 1 after second write: \n");
	for (ii = 0; ii < edges1[0]; ii++) {
	    for (jj = 0; jj < DIM1; jj++) {
		for (kk = 0; kk < DIM2; kk++) {
		fprintf(stderr, " %d ", outdata[ii][jj][kk]);
		}
		fprintf(stderr, "\n");
 		}
 	    fprintf(stderr, "\n");
	}
        fprintf(stderr, "\n\n");

	istat = SDendaccess(sds1_id);

        sds_idx = SDnametoindex (sd_id, "Dataset 2");
        sds2_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds2_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d\n", name1, rank, dimsizes[0], dimsizes[1]);

	/* write the second time to Dataset 2 1x3 */
	start_w[0] = 3;
	edges1[0] = 1;  /* writing second slab */
	istat = SDwritedata(sds2_id, start_w, NULL, edges1, (VOIDP)data2_slab2); 
	CHECKFAIL(istat, " SDwritedata Dataset 2 first time failed \n");

	HDmemset(&outdata2, 0, sizeof(outdata2));
	edges1[0] = 4; /* already wrote 2 slabs and skipped 2 */
        istat = SDreaddata(sds2_id, start, NULL, edges1, (VOIDP) outdata2);
	CHECKFAIL(istat, " SDreaddata Dataset 2 failed \n");
        fprintf(stderr, "\nData read from Dataset 2 after second write: \n");
	for (ii = 0; ii < edges1[0]; ii++) {
	    for (jj = 0; jj < DIM1; jj++) {
		fprintf(stderr, " %d ", outdata2[ii][jj]);
		}
 	    fprintf(stderr, "\n");
	}
        fprintf(stderr, "\n\n");

	istat = SDendaccess(sds2_id);

	/* Display dimension sizes before closing the second time */
        sds_idx = SDnametoindex (sd_id, "Dataset 2");
        sds2_id = SDselect (sd_id, sds_idx);

	istat = SDgetinfo(sds2_id, name1, &rank, dimsizes, &dtype, &nattrs);

 fprintf(stderr, "dataset %s has %d dimensions: %d,%d\n", name1, rank, dimsizes[0], dimsizes[1]);

	istat = SDendaccess(sds2_id);

	/* close the file */
	istat = SDend(sd_id);

	/* Reopen the file for writing. */
	sd_id = SDstart("hdffile_unlim.hdf", DFACC_RDWR);
	if (sd_id == FAIL) {
		fprintf(stderr,"SDstart failed.\n");
		HEprint(stdout,0);
		exit (-1);
	}
	else
		fprintf(stderr,"\n... hdffile_unlim.hdf opened second time\n");

/***** write third and fourth 2x3x2 slabs at once *****/

	/* Fill the array with values for two slab. */
	for (ii = 0; ii < 2; ii++)
	    for (jj = 0; jj < DIM1; jj++)
		for (kk = 0; kk < DIM2; kk++) {
		    data1_slab3_4[ii][jj][kk] =  kk + jj + (ii+3)*100;
		    }

        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds1_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, rank, dimsizes[0], dimsizes[1], dimsizes[2]);

	/* Define the location, pattern, and size of the data set */
	start_w[0] = 4;

	edges1[0] = 2;
	edges1[1] = DIM1;
	edges1[2] = DIM2;
        
	/* write the third time to Dataset 1 */
	istat = SDwritedata(sds1_id, start_w, NULL, edges1, (VOIDP)data1_slab3_4); 
	CHECKFAIL(istat, " SDwritedata Dataset 1 first time failed \n");
 fprintf(stderr, "\n**** After writing third and fourth slabs ****\n");

	HDmemset(&outdata, 0, sizeof(outdata));
	edges1[0] = 6; /* already wrote 6 slabs */
        istat = SDreaddata(sds1_id, start, NULL, edges1, (VOIDP) outdata);
	CHECKFAIL(istat, " SDreaddata failed \n");
        fprintf(stderr, "\nData read from Dataset 1 after third write: \n");
	for (ii = 0; ii < edges1[0]; ii++) {
	    for (jj = 0; jj < DIM1; jj++) {
		for (kk = 0; kk < DIM2; kk++) {
		fprintf(stderr, " %d ", outdata[ii][jj][kk]);
		}
		fprintf(stderr, "\n");
 		}
 	    fprintf(stderr, "\n");
	}
        fprintf(stderr, "\n\n");

	istat = SDendaccess(sds1_id);

        sds_idx = SDnametoindex (sd_id, "Dataset 1");
        sds1_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds1_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d,%d\n", name1, rank, dimsizes[0], dimsizes[1], dimsizes[2]);
	istat = SDendaccess(sds1_id);
	if (istat == -1) fprintf(stderr," SDendaccess failed \n");

        sds_idx = SDnametoindex (sd_id, "Dataset 2");
        sds2_id = SDselect (sd_id, sds_idx);
	istat = SDgetinfo(sds2_id, name1, &rank, dimsizes, &dtype, &nattrs);
 fprintf(stderr, "dataset %s has %d dimensions: %d,%d\n", name1, rank, dimsizes[0], dimsizes[1]);
	istat = SDendaccess(sds2_id);

        istat = SDend(sd_id);
}