/*
 * This file is part of tela the Tensor Language.
 * Copyright (c) 1994-1996 Pekka Janhunen
 */

/*
	CDF.ct
	Notice that CDF is not the same as netCDF, which is already supported by HDF library!
	Preprocess with ctpp.
	C-tela code is C++ equipped with []=f() style function definition.

	---
	This file reads both rVars and zVars but may still have some problems
	with (probably) degenerate cases.
*/

#include "cdf.h"
#include <ctype.h>

#define DEBUG 0

static int StatusHandler(CDFstatus status, int linenumber)
// Return -1 in case of fatal error, 0 otherwise
{
	char message[CDF_ERRTEXT_LEN+1];

	if (status < CDF_WARN) {
		CDFerror(status,message);
		cout << "CDF error " << status << " at line " << linenumber << ":\n  " << message << "\n";
		return -1;
	} else if (status < CDF_OK) {
		cout << "CDF warning " << status << " at line " << linenumber << ":\n  " << message << "\n";
	} else if (status > CDF_OK) {
		CDFerror (status, message);
		cout << "CDF informative message " << status << " at line " << linenumber << ":\n  " << message << "\n";
	}
	return 0;
}

static void MakeValidIdentifier(char *s)
// Replace bad characters in s by underscore (dollar if first character is bad).
// Valid identifier starts with a letter or dollar sign, the rest of the chars
// may contain also digits or underscores. (See lexer source file "d.l").
// Before all, however, the function deletes trailing blanks.
{
	int L = strlen(s) - 1;
	while (s[L] == ' ') {s[L] = '\0'; L--;}
	if (s[0] && !isalpha(s[0]) && s[0]!='$') s[0] = '$';
	int i;
	for (i=1; s[i]; i++)
		if (!isalpha(s[i]) && !isdigit(s[i]) && s[i]!='$' && s[i]!='_')
			s[i] = '_';
}

#define chk(f) status=f; if (status != CDF_OK) {retval=StatusHandler(status,__LINE__); if (retval < 0) return -2;}
#define CDFlong long

static char *DataTypeString(long datatype)
{
	char *s;
	switch (datatype) {
	case CDF_INT1: s = "int1"; break;
	case CDF_INT2: s = "int2"; break;
	case CDF_INT4: s = "int4"; break;
	case CDF_UINT1: s = "uint1"; break;
	case CDF_UINT2: s = "uint2"; break;
	case CDF_UINT4: s = "uint4"; break;
	case CDF_REAL4: s = "real4"; break;
	case CDF_REAL8: s = "real8"; break;
	case CDF_EPOCH: s = "epoch"; break;
	case CDF_BYTE: s = "byte"; break;
	case CDF_FLOAT: s = "float"; break;
	case CDF_DOUBLE: s = "double"; break;
	case CDF_CHAR: s = "char"; break;
	case CDF_UCHAR: s = "uchar"; break;
	default: s = "<UNKNOWN DATA TYPE>";
	}
	return s;
}

// Crash: sample2 cacsst1 d103a1x



static int ReadOneCDF(CDFlong v)
{
	CDFlong rank, dims[CDF_MAX_DIMS], encoding, majority, maxrec, Nvars;
	CDFstatus status;
	int retval, returnvalue=0, d;
	
	char varname[1000];
	CDFlong datatype, Nelems, recvariance, dimvariances[CDF_MAX_DIMS];

	chk(CDFlib(SELECT_,zVAR_,v,NULL_));
	
	chk(CDFlib(GET_,zVAR_NUMDIMS_,&rank,NULL_));
	if (rank > MAXRANK-1) return -4;
	chk(CDFlib(GET_,zVAR_DIMSIZES_,dims,NULL_));
	chk(CDFlib(GET_,zVAR_DATATYPE_,&datatype,NULL_));
	chk(CDFlib(GET_,zVAR_NAME_,varname,NULL_));
	chk(CDFlib(GET_,zVAR_MAXREC_,&maxrec,NULL_));
	chk(CDFlib(GET_,zVAR_NUMELEMS_,&Nelems,NULL_));
	chk(CDFlib(GET_,zVAR_RECVARY_,&recvariance,NULL_));
	if (rank == 0) {rank = 1; dims[0] = 1;}
		
#if DEBUG
	cout << "rank=" << rank << ", dims=";
	for (d=0; d<rank; d++) cout << " " << dims[d];
	cout << "- Variable '" << varname << "', Nelems="
		 << Nelems << ", datatype=" << DataTypeString(datatype) << ", recvariance=" << recvariance << "\n" << flush;
#endif

	long itemsize;
	switch (datatype) {
	case CDF_INT1: case CDF_BYTE: case CDF_CHAR: case CDF_UCHAR: itemsize = sizeof(char); break;
	case CDF_INT2: itemsize = sizeof(short int); break;
	case CDF_INT4: itemsize = sizeof(int); break;
	case CDF_UINT1: itemsize = sizeof(unsigned char); break;
	case CDF_UINT2: itemsize = sizeof(unsigned short int); break;
	case CDF_UINT4: itemsize = sizeof(unsigned int); break;
	case CDF_REAL4: case CDF_FLOAT: itemsize = sizeof(float); break;
	case CDF_REAL8: case CDF_DOUBLE: case CDF_EPOCH: itemsize = sizeof(double); break;
	default: itemsize = -1;
		
	}

	if (itemsize == -1) {
		cout << "---- Warning: Unsupported data type\n" << flush;
		return -10;
	}
	long Nrecs = (recvariance == VARY) ? maxrec : 1;
	long recStart = 0;
	long recCount = (recvariance == VARY) ? maxrec : 1;
	long recInterval = 1;
	long *dimStart = new long [3*rank];
	long *dimEnd = dimStart + rank;
	long *dimIntervals = dimStart + 2*rank;
	long length = 1;
	Tint objdims[MAXRANK];
#if DEBUG
	cout << "--- Nrecs=" << Nrecs << "\n" << flush;
#endif
		
	if (Nrecs > 1) objdims[0] = Nrecs;
#if DEBUG
	cout << "--- rank=" << rank << "\n" << flush;
#endif
	dimStart[0] = dimEnd[0] = 0; dimIntervals[0] = 1;
	if (rank != 1 || dims[0] != 1 || Nelems > 1 || Nrecs > 1 || (rank==1 && dims[0]==1 && Nelems==1 && Nrecs==1)) {
		for (d=0; d<rank; d++) {
			dimStart[d] = 0;
			dimEnd[d] = dims[d];
			objdims[d + (Nrecs > 1 ? 1 : 0)] = Tint(dims[d]);
			length*= dims[d];
			dimIntervals[d] = 1;
		}
	}
	if (Nelems > 1) {
		objdims[d + (Nrecs > 1 ? 1 : 0)] = Nelems;
		length*= Nelems;
	}
	TDimPack odims(objdims,
				   int((rank==1 && dims[0]==1 && (Nelems!=1 || Nrecs!=1) ? 0 : rank)
					   + (Nrecs > 1 ? 1 : 0)
					   + (Nelems > 1 ? 1 : 0) ));

#if DEBUG
	cout << "odims = " << odims << "\n";
	cout << "itemsize=" << itemsize << ", length=" << length << ", recCount=" << recCount << "\n";
	cout << "Allocating buffer of size=" << itemsize*length*recCount << "\n" << flush;
#endif
	char *buffer = new char [itemsize*length*recCount];

	chk(CDFlib(SELECT_,
			   zVAR_,v,
			   zVAR_RECNUMBER_, recStart,
			   zVAR_RECCOUNT_, recCount,
			   zVAR_RECINTERVAL_, recInterval,
			   zVAR_DIMINDICES_, dimStart,
			   zVAR_DIMCOUNTS_, dimEnd,
			   zVAR_DIMINTERVALS_, dimIntervals,
			   GET_, zVAR_HYPERDATA_, buffer,
			   NULL_));

	Tobject obj;
	Treal *rptr; Tint *iptr;
	switch (datatype) {
		// All integer types -> Tint, except UINT4 which maps to Treal
		// because UINT4 will not fit into Tint on 32-bit machines!
	case CDF_INT1:
	case CDF_BYTE:
	case CDF_INT2:
	case CDF_INT4:
	case CDF_UINT1:
	case CDF_UINT2:
	case CDF_CHAR:
	case CDF_UCHAR:
		obj.izeros(odims); iptr = obj.IntPtr(); break;
	case CDF_REAL4:
	case CDF_REAL8:
	case CDF_FLOAT:
	case CDF_DOUBLE:
	case CDF_EPOCH:
	case CDF_UINT4:
		obj.rzeros(odims); rptr = obj.RealPtr(); break;
	default:
		cout << "*** Unsupported data type!\n" << flush;
		delete [] buffer;
		delete [] dimStart;
		return -10;
	}
	
	Tint i, N = obj.length()  /*Tint(length*Nrecs)*/;
	switch (datatype) {
	case CDF_INT1:
	case CDF_BYTE:
	case CDF_CHAR:
		for (i=0; i<N; i++) iptr[i] = ((char *)buffer)[i];
		if (datatype == CDF_CHAR) obj.SetStringFlag();
		break;
	case CDF_INT2:
		for (i=0; i<N; i++) iptr[i] = ((short int*)buffer)[i];
		if (sizeof(short int) != 2) returnvalue = 1;
		break;
	case CDF_INT4:
		for (i=0; i<N; i++) iptr[i] = ((int*)buffer)[i];
		if (sizeof(int) != 4) returnvalue = 2;
		break;
	case CDF_UINT1:
	case CDF_UCHAR:
		for (i=0; i<N; i++) iptr[i] = ((unsigned char*)buffer)[i];
		if (datatype == CDF_UCHAR) obj.SetStringFlag();
		break;
	case CDF_UINT2:
		for (i=0; i<N; i++) iptr[i] = ((unsigned short int*)buffer)[i];
		if (sizeof(unsigned short int) != 2) returnvalue = 1;
		break;
	case CDF_UINT4:
		for (i=0; i<N; i++) rptr[i] = ((unsigned int*)buffer)[i];
		if (sizeof(unsigned int) != 4) returnvalue = 2;
		break;
	case CDF_REAL4:
	case CDF_FLOAT:
		for (i=0; i<N; i++) rptr[i] = ((float *)buffer)[i];
		if (sizeof(float) != 4) returnvalue = 3;
		break;
	case CDF_REAL8:
	case CDF_DOUBLE:
	case CDF_EPOCH:
		for (i=0; i<N; i++) rptr[i] = ((double *)buffer)[i];
		if (sizeof(double) != 8) returnvalue = 4;
		break;
	default: cout << "Unsupported type!!!!!\n" << flush;
	}
	
	MakeValidIdentifier(varname);
#if DEBUG
	cout << "After MakeValidIdentifier, varname='" << varname << "'\n";
#endif
	if (strlen(varname) > 0) {
		Tsymbol *symptr = theHT.add((Tchar*)varname);
		symptr->value()->bitwiseMoveFrom(obj);
	} else {
		cout << "- Empty variable name in CDF file - ignoring\n";
		returnvalue = 5;
	}
	
	delete [] buffer;
	delete [] dimStart;
#if DEBUG
	cout << flush;
#endif
	return returnvalue;
}



[] = import_CDF(fn...)
/* import_CDF("fyle") reads the variables in CDF file "fyle.cdf"
   into Tela workspace. Notice that you must leave out the suffix .cdf!

   import_CDF("fyle","var1","var2",...) reads only variables with
   specified names.
   
   Limitations:
   - Only rVariables are currently recognized
   - The row major/col major flag is ignored
   Error codes:
   1: Loading INT2/UINT2 CDF data, but this Tela kernel has sizeof(short int) != 2
   2: Loading INT4/UINT4 CDF data, but this Tela kernel has sizeof(int) != 4
   3: Loading FLOAT CDF data, but this Tela kernel has sizeof(float) != 4
   4: Loading DOUBLE CDF data, but this Tela kernel has sizeof(double) != 8
   5: Empty variable name in CDF file - ignoring
   -1: First argument not a string
   -2: Fatal CDF library error
   -3: Unsupported data type
   -4: Too high rank in CDF file for this Tela kernel
   -5: Function called with more than one argument and the rest are not all strings
   -10: Internal error
*/
{
	if (!fn.IsString()) return -1;
	Tstring FN = fn;
	char *f = (char*)FN;
	int d,retval,returnvalue=0;
	CDFid id;
	CDFstatus status;
	CDFlong v,rank, dims[CDF_MAX_DIMS], encoding, majority, maxrec, Nvars;
	
	chk(CDFopen(f,&id));
	chk(CDFlib(SELECT_,
			   CDF_READONLY_MODE_,READONLYon,
			   CDF_zMODE_,zMODEon2,
			   CDF_,id,
			   NULL_));
	if (Nargin > 1) {

		int a;
		for (a=1; a<Nargin; a++) {
			if (!argin[a]->IsString()) return -5;
			Tstring S = *argin[a];
#if DEBUG
			cout << "S = " << S << "\n" << flush;
#endif
			chk(CDFlib(GET_,zVAR_NUMBER_,(char*)S,&v,NULL_));
#if DEBUG
			cout << "- Variable name '" << (char*)S << "' has number " << v << "\n" << flush;
#endif
			returnvalue = ReadOneCDF(v);
			if (returnvalue < 0) {CDFclose(id); return returnvalue;}
		}
		
	} else {
		chk(CDFlib(GET_,CDF_NUMzVARS_,&Nvars,NULL_));
#if DEBUG
		cout << "Nvars = " << Nvars << "\n";
#endif
		for (v=0; v<Nvars; v++) {
			chk(CDFlib(SELECT_,zVAR_,v,NULL_));
			returnvalue = ReadOneCDF(v);
			if (returnvalue < 0) {CDFclose(id); return returnvalue;}
		}
	}
	
	chk(CDFclose(id));
	
	return returnvalue;
}















#if 0

static void foo()
{
	if (!fn.IsString()) return -1;
	Tstring FN = fn;
	char *f = (char*)FN;
	int d,v,retval,returnvalue=0;
	CDFid id;
	CDFstatus status;
	CDFlong rank, dims[CDF_MAX_DIMS], encoding, majority, maxrec, Nvars, Nattrs;
	
	chk(CDFopen(f,&id));
	chk(CDFinquire(id,&rank,dims,&encoding,&majority,&maxrec,&Nvars,&Nattrs));

#if DEBUG
	cout << "CDFinquire: rank=" << rank << ", dims=";
	for (d=0; d<rank; d++) cout << " " << dims[d];
	cout << ", encoding=" << encoding;
	cout << ", majority=" << majority << ", maxrec=" << maxrec << ", Nvars=" << Nvars << ", Nattrs=" << Nattrs << "\n";
#endif

	if (rank > MAXRANK-1) return -4;

	for (v=0; v<Nvars; v++) {

		char varname[1000];
		CDFlong datatype, Nelems, recvariance, dimvariances[CDF_MAX_DIMS];
		chk(CDFvarInquire(id,v,varname,&datatype,&Nelems,&recvariance,dimvariances));
#if DEBUG
		cout << "- Variable " << v << ": name='" << varname << "', Nelems="
			 << Nelems << ", datatype=" << DataTypeString(datatype) << ", recvariance=" << recvariance << "\n" << flush;
#endif

		long itemsize;
		switch (datatype) {
		case CDF_INT1: case CDF_BYTE: case CDF_CHAR: case CDF_UCHAR: itemsize = sizeof(char); break;
		case CDF_INT2: itemsize = sizeof(short int); break;
		case CDF_INT4: itemsize = sizeof(int); break;
		case CDF_UINT1: itemsize = sizeof(unsigned char); break;
		case CDF_UINT2: itemsize = sizeof(unsigned short int); break;
		case CDF_UINT4: itemsize = sizeof(unsigned int); break;
		case CDF_REAL4: case CDF_FLOAT: itemsize = sizeof(float); break;
		case CDF_REAL8: case CDF_DOUBLE: case CDF_EPOCH: itemsize = sizeof(double); break;
//		case CDF_EPOCH: itemsize = 1; Nelems = EPOCH_STRING_LEN; break;
		default: itemsize = -1;
			
		}

		if (itemsize == -1) {
			cout << "---- Warning: Unsupported data type\n" << flush;
			continue;
		}
		long Nrecs = (recvariance == VARY) ? maxrec : 1;
		long recStart = 0;
		long recCount = (recvariance == VARY) ? maxrec : 1;
		long recInterval = 1;
		long *dimStart = new long [3*rank];
		long *dimEnd = dimStart + rank;
		long *dimIntervals = dimStart + 2*rank;
		long length = 1;
		Tint objdims[MAXRANK];
#if DEBUG
		cout << "--- Nrecs=" << Nrecs << "\n" << flush;
#endif
		
		if (Nrecs > 1) objdims[0] = Nrecs;
#if DEBUG
		cout << "--- rank=" << rank << "\n" << flush;
#endif
		for (d=0; d<rank; d++) {
			dimStart[d] = 0;
			dimEnd[d] = dims[d];
			objdims[d + (Nrecs > 1 ? 1 : 0)] = Tint(dims[d]);
			length*= dims[d];
			dimIntervals[d] = 1;
		}
		if (Nelems > 1) {
			objdims[d + (Nrecs > 1 ? 1 : 0)] = Nelems;
			length*= Nelems;
		}
		TDimPack odims(objdims,int(rank + (Nrecs > 1 ? 1 : 0) + (Nelems > 1 ? 1 : 0)));

#if DEBUG
		cout << "odims = " << odims << "\n";
		cout << "itemsize=" << itemsize << ", length=" << length << ", recCount=" << recCount << "\n";
		cout << "Allocating buffer of size=" << itemsize*length*recCount << "\n" << flush;
#endif
		char *buffer = new char [itemsize*length*recCount];

		chk(CDFvarHyperGet(id,v,recStart,recCount,recInterval,dimStart,dimEnd,dimIntervals,buffer));

		Tobject obj;
		Treal *rptr; Tint *iptr;
		switch (datatype) {
			// All integer types -> Tint, except UINT4 which maps to Treal
			// because UINT4 will not fit into Tint on 32-bit machines!
		case CDF_INT1:
		case CDF_BYTE:
		case CDF_INT2:
		case CDF_INT4:
		case CDF_UINT1:
		case CDF_UINT2:
		case CDF_CHAR:
		case CDF_UCHAR:
			obj.izeros(odims); iptr = obj.IntPtr(); break;
		case CDF_REAL4:
		case CDF_REAL8:
		case CDF_FLOAT:
		case CDF_DOUBLE:
		case CDF_EPOCH:
		case CDF_UINT4:
			obj.rzeros(odims); rptr = obj.RealPtr(); break;
		default:
			cout << "*** Unsupported data type!\n" << flush;
			delete [] buffer;
			delete [] dimStart;
			continue;
		}

		Tint i, N = Tint(length*Nrecs);
		switch (datatype) {
		case CDF_INT1:
		case CDF_BYTE:
		case CDF_CHAR:
			for (i=0; i<N; i++) iptr[i] = ((char *)buffer)[i];
			if (datatype == CDF_CHAR) obj.SetStringFlag();
			break;
		case CDF_INT2:
			for (i=0; i<N; i++) iptr[i] = ((short int*)buffer)[i];
			if (sizeof(short int) != 2) returnvalue = 1;
			break;
		case CDF_INT4:
			for (i=0; i<N; i++) iptr[i] = ((int*)buffer)[i];
			break;
		case CDF_UINT1:
		case CDF_UCHAR:
			for (i=0; i<N; i++) iptr[i] = ((unsigned char*)buffer)[i];
			if (datatype == CDF_UCHAR) obj.SetStringFlag();
			break;
		case CDF_UINT2:
			for (i=0; i<N; i++) iptr[i] = ((unsigned short int*)buffer)[i];
			if (sizeof(unsigned short int) != 2) returnvalue = 1;
			break;
		case CDF_UINT4:
			for (i=0; i<N; i++) iptr[i] = ((unsigned int*)buffer)[i];
			break;
		case CDF_REAL4:
		case CDF_FLOAT:
			for (i=0; i<N; i++) rptr[i] = ((float *)buffer)[i];
			break;
		case CDF_REAL8:
		case CDF_DOUBLE:
		case CDF_EPOCH:
			for (i=0; i<N; i++) rptr[i] = ((double *)buffer)[i];
			break;
		default: cout << "Unsupported type!!!!!\n" << flush;
		}

		MakeValidIdentifier(varname);
		Tsymbol *symptr = theHT.add((Tchar*)varname);
		symptr->value()->bitwiseMoveFrom(obj);

		delete [] buffer;
		delete [] dimStart;
#if DEBUG
		cout << flush;
#endif
	}
	
	chk(CDFclose(id));
	
	return returnvalue;

}

#endif