Bug fixes in SPICE rotations

Provides various bug fixes associated with handling non-cartesian
vector data including:
* vector values not mapped to canonical directions on write
* output rotation data shape is union of time variable and input
  variable
* Extra coordinates not mapped to a plot axis are allowed

das2/dataset_hdr3.c

@@ -73,7 +73,8 @@ typedef struct serial_xml_context {
 	char valStorage[_VAL_STOREAGE_SZ];
 	ubyte varCompSys;
 	ubyte varCompDirs;
-	int nVarComps;
+	int nVarComps;  /* only Non-zero if the item is a vector.  
+	                   A scalar is *not* a 1-component vector! */
 
 	char varCompLbl[_VAL_COMP_LBL_SZ]; /* HACK ALERT: Temporary hack for dastelem output */
 
@@ -344,10 +345,12 @@ static void _serial_onOpenDim(
 		if(strcmp(psAttr[i],"physDim")==0)     sPhysDim = psAttr[i+1];
 		else if(strcmp(psAttr[i],"name")==0)   sName    = psAttr[i+1];
 		else if(strcmp(psAttr[i],"frame")==0)  sFrame   = psAttr[i+1];
-		else if((strcmp(psAttr[i],"axis")==0) &&(psAttr[i+1][0] != '\0')) 
-			strncpy(sAxis, psAttr[i+1], 47);
-		else if((strcmp(psAttr[i],"annotation")==0) &&(psAttr[i+1][0] != '\0')) 
-			strncpy(sAnnot, psAttr[i+1], 47);
+		else if(strcmp(psAttr[i],"axis")==0){ 
+			if(psAttr[i+1][0] != '\0') strncpy(sAxis, psAttr[i+1], 47);
+		}
+		else if(strcmp(psAttr[i],"annotation")==0){
+			if(psAttr[i+1][0] != '\0') strncpy(sAnnot, psAttr[i+1], 47);
+		} 
 		else
 			daslog_warn_v(
 				"Unknown attribute %s in <%s> for dataset ID %02d", psAttr[i], sDimType, id
@@ -367,13 +370,15 @@ static void _serial_onOpenDim(
 	if(sPhysDim[0] == '\0')
 		sPhysDim = "none";
 
+	/* 
 	if((dt == DASDIM_COORD) && (sAxis[0] == '\0') && (sAnnot[0] == '\0')){
 		pCtx->nDasErr = das_error(DASERR_SERIAL, 
 			"Both \"axis\" and \"annotation\" missing for coordinate dimension %s in dataset ID %d", 
 			sPhysDim, id
 		);
 		return;
 	}
+	*/
 
 	/* We have required items, make the dim */
 	DasDim* pDim = new_DasDim(sPhysDim, sName, dt, DasDs_rank(pCtx->pDs));
@@ -736,7 +741,7 @@ static DasErrCode _serial_makeVarAry(context_t* pCtx, bool bHandleFill)
 
 	if(pCtx->varIntRank > 0){	
 		/* Internal structure due to vectors */
-		if(pCtx->nVarComps > 1){
+		if(pCtx->nVarComps > 0){
 			aShape[nAryRank] = pCtx->nVarComps;
 			++nAryRank;
 		}

das2/dimension.c

@@ -416,7 +416,7 @@ DasErrCode DasVar_encode(DasVar* pVar, const char* sRole, DasBuf* pBuf);
 DasErrCode DasDim_encode(DasDim* pThis, DasBuf* pBuf)
 {
 	DasErrCode nRet = DAS_OKAY;
-	char sAxis[64];
+	char sAxis[64] = {'\0'};
 	const char* sType = "data";
 	if(DasDim_type(pThis) == DASDIM_COORD){
 		sType = "coord";
@@ -436,7 +436,7 @@ DasErrCode DasDim_encode(DasDim* pThis, DasBuf* pBuf)
 		sType, DasDim_dim(pThis), DasDim_id(pThis)
 	);
 
-	if(DasDim_type(pThis) == DASDIM_COORD){
+	if((DasDim_type(pThis) == DASDIM_COORD)&&(sAxis[0] != '\0')){
 		if(DasDim_isPrimeCoord(pThis))
 			DasBuf_printf(pBuf, " axis=\"%s\"", sAxis);
 		else

das2/variable.h

@@ -464,7 +464,7 @@ DAS_API DasVar* new_DasConstant(const char* sId, const das_datum* pDm);
 /** Create a simple linear sequence variable
  * 
  * A simple sequence variable is linear in a single index.  Many measurements
- * happen in a parameter that progresses as a simple linea function of a 
+ * happen in a parameter that progresses as a simple linear function of a 
  * single index.  For example time offest for a single A/D capture from the
  * start of the capture sequence.
  * 
@@ -516,8 +516,10 @@ DAS_API DasVar* new_DasVarSeq(
  *          Don't set this directly if you can avoid it.  Use the 
  *          SCALAR_N and VECTOR_N macros instead.
  * 
- * @param pMap The mapping of external indexs to DasAry indexes.  Not every
- *          external index needs to be mapped.
+ * @param pMap The mapping of external indexes to DasAry indexes.  The offset
+ *          into this array is the external index.  The value is the internal
+ *          index.  Negative values indicate that the is *no* external mapping.
+ *          Not every external index needs to be mapped.
  *
  *          Don't set this directly if you can avoid it.  Use the SCALAR_N
  *          and VECTOR_N macros instead.

das2/vector.c

@@ -267,55 +267,70 @@ DasErrCode das_geovec_values(das_geovec* pThis, double* pValues)
 
 	size_t i = 0;
 
+	/* set default values to handle missing components.  Only matters
+		for non-cartesian systems */
+
+	if((pThis->systype == DAS_VSYS_CYL)||(pThis->systype == DAS_VSYS_SPH)||
+		(pThis->systype == DAS_VSYS_CENTRIC)
+	)
+		pValues[0] = 1.0;
+
+	/* Remap based on dirs, sure wish I had room to save this away */
+	int dirs[3] = {0};
+	for(int i = 0; i < pThis->ncomp; ++i)
+		dirs[i] = (pThis->dirs << i*2)&0x3;	
+
 	switch(pThis->et){
 	case vtByte:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((int8_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((int8_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtUByte:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((ubyte*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((ubyte*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtShort:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((int16_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((int16_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtUShort:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((uint16_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((uint16_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtInt:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((int32_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((int32_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtUInt:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((uint32_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((uint32_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtLong:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((int64_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((int64_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtULong:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((uint64_t*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((uint64_t*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtFloat:
 		for(i = 0; i < pThis->ncomp; ++i)
-			pValues[i] = ((float*)(pThis->comp))[ i ];
+			pValues[dirs[i]] = ((float*)(pThis->comp))[ i ];
 		return DAS_OKAY;
 
 	case vtDouble:
-		memcpy(pValues, pThis->comp, (pThis->ncomp)*sizeof(double));
+		for(i = 0; i < pThis->ncomp; ++i)
+			pValues[dirs[i]] = ((double*)(pThis->comp))[ i ];
 		return DAS_OKAY;
+
 	default: break;
 	}
 

utilities/das3_spice.c

@@ -872,34 +872,55 @@ DasErrCode _addRotation(XCalc* pCalc, const char* sAnonFrame, DasDs* pDsOut)
 		strncpy(pReq->aInFrame, sFrame, DASFRM_NAME_SZ - 1);
 	}
 
-	/* The shape the storage array is just the same as the input, with all
-		unused indexes collapsed.  A three-way compare is used in order to
-		support function variables, which take the same shape as thier 
-		container */
+	/* The shape the storage array is:
+
+	   Shape of the input + shape of the time reference
+	*/
 	ptrdiff_t aVarShape[DASIDX_MAX] = DASIDX_INIT_UNUSED;
 	DasVar_shape(pCalc->pVarIn, aVarShape);
 
+	ptrdiff_t aTimeShape[DASIDX_MAX] = DASIDX_INIT_UNUSED;
+	DasVar_shape(pCalc->pTime, aTimeShape);
+
+	/* Take the union of the shapes */
+	ptrdiff_t aCombined[DASIDX_MAX] =  DASIDX_INIT_UNUSED;
+	das_varindex_merge(DASIDX_MAX - 1, aCombined, aTimeShape);
+	das_varindex_merge(DASIDX_MAX - 1, aCombined, aVarShape);
+
+	/* Now the array shape is all the used items, but watch out
+	   for ragged */
 	int nAryRank = 0;
-	size_t aAryShape[DASIDX_MAX] = {0};
-	int i,j;
 	int nItems = 1;
+	size_t aAryShape[DASIDX_MAX] = {0};
+
+	/* the index into this array is the the overall dataset external
+	   index.  The values represent the array "dimension" */
 	int8_t aVarMap[DASIDX_MAX] = DASIDX_INIT_UNUSED;
-	for(i = 0, j = 0; i < nDsRank; ++i){
-		if(!DasVar_degenerate(pCalc->pVarIn, i)) aVarMap[i] = (int8_t)i;
 
-		/* Even if upstream doesn't use the record index, pTime does, so we do too */
-		if((i > 0) && (aVarMap[i] == DASIDX_UNUSED)) continue;
+	for(int iExtern = 0; iExtern < DASIDX_MAX; ++iExtern){
 
-		aAryShape[j] = (aDsShape[i] == DASIDX_RAGGED) ? 0 : aDsShape[i];
+		// If unused I have no internal array index value... */
+		if(aCombined[iExtern] == DASIDX_UNUSED)
+			continue;
+
+		aVarMap[iExtern] = nAryRank;  // ...otherwise map increasing
+
+		aAryShape[nAryRank] = 
+			(aCombined[iExtern] == DASIDX_RAGGED) ? 0 : aCombined[iExtern];
+
+		assert(aAryShape[nAryRank] >= 0);
 
-		/* The items per record are not affected by the record index */
-		if(i > 0) nItems *= aAryShape[j];
+		/* Items per record matters to the codec, so save that */
+		if(iExtern > 0) 
+			nItems *= aAryShape[nAryRank];
 
 		++nAryRank;
 	}
+
 	aAryShape[nAryRank] = 3;  /* Add 1 internal dimension */
 	nItems *= 3;
 	++nAryRank;
+
 	if(nItems < 0) nItems = -1; /* if any dim ragged, just make codec ragged */
 
 	const DasDim* pDimIn = (const DasDim*)DasDesc_parent((const DasDesc*)(pCalc->pVarIn));
@@ -1400,6 +1421,8 @@ DasErrCode _writeRotation(DasDs* pDsIn, XCalc* pCalc, double rTimeShift)
 		das_geovec_values(pVecIn, aRecIn);
 
 		if(pVecIn->systype != DAS_VSYS_CART){   /* convert non-cart input coords */
+			memcpy(aTmp, aRecIn, sizeof(SpiceDouble)*3);
+
 			switch(pVecIn->systype){
 
 			case DAS_VSYS_CYL:          /* Assume degrees, but need to check */