heap.c

Stack-based sequential decoder for M-QAM modulated MIMO-type problems, i.e., of fixed tree depth.

#ifdef DEBUG
      printf("SwapRight "); 
#endif

    }

    /* link appropriate child (TempNode) to parent or root */
    if (pNewNode->pPar) { 
      if (pNewNode->pPar->pLCh == pNewNode) { /* NewNode was left child */
        pNewNode->pPar->pLCh = pTempNode;
      } else {                                /* NewNode was right child */
        pNewNode->pPar->pRCh = pTempNode;
      }
      pTempNode->pPar = pNewNode->pPar;
    } else {
      pH->pRoot = pTempNode;
      pTempNode->pPar = NULL;
    }

    /* Finish linking appropriate child as parent */
    pNewNode->pPar = pTempNode; 

    /* link left grandchild as left child */
    pNewNode->pLCh = pTempLCh;
    if (pTempLCh) { pTempLCh->pPar = pNewNode; }

    /* link right grandchild as right child */
    pNewNode->pRCh = pTempRCh;
    if (pTempRCh) { pTempRCh->pPar = pNewNode; }

    /* recompute loop condition */
    swapLCh = pNewNode->pLCh && 
	        ( pNewNode->pLCh->key < pNewNode->key ||
	          ( pNewNode->pLCh->key == pNewNode->key && 
	            pNewNode->pLCh->dim < pNewNode->dim)  );
    swapRCh = pNewNode->pRCh &&
	        ( pNewNode->pRCh->key < pNewNode->key ||
	          ( pNewNode->pRCh->key == pNewNode->key && 
	            pNewNode->pRCh->dim < pNewNode->dim)  );
    if (swapLCh && swapRCh) {
      swapLCh = ( pNewNode->pLCh->key < pNewNode->pRCh->key ||
	          ( pNewNode->pLCh->key == pNewNode->pRCh->key && 
	            pNewNode->pLCh->dim <= pNewNode->pRCh->dim)  );
      swapRCh = !swapLCh;
    } 
  }  /* while */

#ifdef DEBUG
  printf("\n"); 
#endif

  return nSwaps;
}

void hNodeFree(struct hNode *pNode) {
  if (pNode->pRCh) {
    hNodeFree(pNode->pRCh);
  }
  if (pNode->pLCh) {
    hNodeFree(pNode->pLCh);
  }
  myfree(pNode->pY);
  myfree(pNode);
}

void hNodeLFree(struct hNodeL *pNodel) {
  if (pNodel->pRCh) {
    hNodeLFree((struct hNodeL *)pNodel->pRCh);
  }
  if (pNodel->pLCh) {
    hNodeLFree((struct hNodeL *)pNodel->pLCh);
  }
  myfree(pNodel->pY);
  if (pNodel->dim>0) { myfree(pNodel->pA); }
  myfree(pNodel);
}

void heapClear(struct heap *pH) {
  if (pH->pRoot) {
    switch (pH->type) {
      case BASIC:
      case LATTICE:
        hNodeFree(pH->pRoot);
        break;
      case LAST:
        hNodeLFree((struct hNodeL *)pH->pRoot);
        break;
      default:
        break;
    }
    pH->pRoot  = NULL;
    pH->nNodes = 0;
  }
}

void heapDelete(struct heap *pH) {
  if (pH->pRoot) {
    switch (pH->type) {
      case BASIC:
      case LATTICE:
        hNodeFree(pH->pRoot);
        break;
      case LAST:
        hNodeLFree((struct hNodeL *)pH->pRoot);
        break;
      default:
        break;
    }
  }
  myfree(pH);
}

void hNodePrintContents(struct hNode *pNode, int dim, int spdim) {
  int ii, jj;

  if (pNode && 
      (pNode->key || pNode->pLCh || pNode->pY || pNode->dim)) {
    for (ii = 0; ii < dim; ii++) {                  /* tabbing */
      printf("     ");
    }
    printf("%5.2f (%i) [", pNode->key, pNode->dim); /* key and dim */
    for (jj = 0; jj < pNode->dim; jj++) {
      printf(" %5.2f ", pNode->pY[jj]);               /* y */
    }		
    printf("] [");
    for (jj = 0; jj < spdim-pNode->dim; jj++) {
      printf(" %5.2f ", pNode->pZ[jj]);               /* z */
    }		
    printf("]\n");

    if (pNode->pLCh) {
      hNodePrintContents(pNode->pLCh, dim+1, spdim);
    } else {
      printf("\n");
    }
    
    if (pNode->pRCh) {
      hNodePrintContents(pNode->pRCh, dim+1, spdim);
    }

  } else {
    printf("[Empty hNode]\n");
  }
}

void heapPrintContents(struct heap *pH) {
  if (pH) {
    if (pH->nNodes > 0) {
      hNodePrintContents(pH->pRoot, 0, pH->spdim); 
    } else {
      printf("[Empty heap]\n");
    } 
  } else {
    printf("[Null heap]\n");
  }
}

int heapGetSize(struct heap *pH) {
  return pH->nNodes;
}

int heapGetMaxSize(struct heap *pH) {
  return pH->maxNodes;
}

int heapGetDim(struct heap *pH) {
  return pH->spdim;
}

void heapPrintUsage() {
  myprintf("\nUsage:\n\n");
  myprintf("       myHeap = heap('Init', problemDimension [, heapSize], type);\n\n");
  myprintf("       nSwaps = heap('Insert', myHeap,\n");
  myprintf("                     biasedNodeWeight,\n");
  myprintf("                     parentResidualTarget,\n");
  myprintf("                     accumulatedSymbolDecisions,\n");
  myprintf("                     residualDimension\n");
  myprintf("                     [, childOrdinal,\n");
  myprintf("                     unbiasedParentNodeWeight,\n");
  myprintf("                     unconstrainedTarget,\n");
  myprintf("                     previousSiblingDecision\n");
  myprintf("                     [, parentShapingSphereOffset,\n");
  myprintf("                     shapingSphereSquaredRadius,\n");
  myprintf("                     siblingRangeLowerBound,\n");
  myprintf("                     siblingRangeUpperBound,\n");
  myprintf("                     unconstrainedOffset]]);\n\n");
  myprintf("       nSwaps = heap('ReplaceMin', myHeap,\n");
  myprintf("                     bNW, pRT, aSD, rD\n");
  myprintf("                     [, cO, uPNW, uT, pSD\n");
  myprintf("                     [, pSSO, sSSR, sRLB, sRUB, uO]]);\n\n");
  myprintf("       [bNW, pRT, aSD, rD [, cO, uPNW, uT, pSD [, pSSO, sSSR, sRLB, sRUB, uO]]]\n");
  myprintf("              = heap('GetMin', myHeap);\n\n");
  myprintf("       nSwaps = heap('DelMin', myHeap);\n\n");
  myprintf("       nNodes = heap('GetSize', myHeap);\n");
  myprintf("       T      = heap('GetMaxSize', myHeap);\n");
  myprintf("       m      = heap('GetDim', myHeap);\n\n");
  myprintf("                heap('Clear', myHeap);\n");
  myprintf("                heap('Delete', myHeap, 'myHeap');\n");
  myprintf("                heap('Print', myHeap);\n\n");
  myprintf("Notes:\n\n       Input vectors may be columns or rows; output vectors are always columns.\n");
  myprintf("       Operations heap('Insert', ...) and heap('ReplaceMin', ...) return -1 when\n");
  myprintf("       out of memory.  Please send bug reports to [e-mail].\n\n");
  myerrmsg("Invalid heap operation, please try again.\n\n");
}

#ifdef COMPILEFORMATLAB
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
  struct heap   *pH;
  struct hNode  *pNode;
  struct hNodeI *pNodei;
  struct hNodeL *pNodel;
  char           Op[20];
  double         key, *pY, *pZ, pweight, utarget, delta, *pA, sqradius, uoffset;
  int            spdim, T, dim, cnum, ns, nn, lbound, ubound, type, nY, nZ;
  mxArray       *tempArray;

  if (nrhs < 1) { 
    heapPrintUsage(); 
    return;
  }

  /* Get the operation */
  mxGetString(prhs[0],Op,20);

  if (mystrncasecmp(Op,"Init",3) == 0) {
    if (nlhs != 1) {
      type = -1;
    } else if (nrhs == 3) {
      type = (int)mxGetScalar(prhs[2]);
    } else if (nrhs == 4) {
      T    = (int)mxGetScalar(prhs[2]);
      type = (int)mxGetScalar(prhs[3]);
    } else {
      type = -1;
    }

    switch (type) {
      case BASIC:
      case LATTICE:
      case LAST:
        spdim = (int)mxGetScalar(prhs[1]);
        if (nrhs == 3) { pH = heapInit(spdim,type); }
        if (nrhs == 4) { pH = heapInitT(spdim,T,type); }
        break;

      default:
        myprintf("\nUsage: myHeap = heap('Init', problemDimension [, heapSize], type);\n\n");
        myprintf("       Initializes a heap object to handle a decoding problem in pD\n");      
        myprintf("       real dimensions using the OTD algorithm.  Note that myHeap\n");
        myprintf("       will appear in Matlab as an empty matrix.  If the heapSize is\n");      
        myprintf("       not specified, myHeap will have a default maxSize of 32768.\n");      
        myprintf("       Set type == 0 for a lattice decoding problem or type == 1 for\n");      
        myprintf("       a Lattice Space-Time (LAST) decoding problem.\n\n");      
        myerrmsg("Invalid heap operation, please try again.\n\n");
        break;
    }

    plhs[0] = mxCreateDoubleMatrix(0,0,mxREAL);
    mxFree(mxGetPr(plhs[0]));
    mxSetPr(plhs[0],(double *)pH);

  } else if (mystrncasecmp(Op,"Insert",3) == 0) {
    if ((nrhs < 6) || (nlhs != 1)) {
      type = -1;
    } else {
      pH = (struct heap *)mxGetPr(prhs[1]);
      type = pH->type;
      key  = (double)mxGetScalar(prhs[2]);
      pY   = mxGetPr(prhs[3]);
      nY   = mxGetN(prhs[3])*mxGetM(prhs[3]);
      pZ   = mxGetPr(prhs[4]);
      nZ   = mxGetN(prhs[4])*mxGetM(prhs[4]);
      dim  = (int)mxGetScalar(prhs[5]);

      if ((nY != dim) || (nY+nZ != pH->spdim)) {
        myprintf("\nError: Either dim(parentResidualTarget) ~= residualDimension or\n");
        myprintf("       dim(pRT)+dim(accumulatedSymbolDecisions) ~= problemDimension.\n");
        myprintf("       All vectors should be row vectors.\n\n");
        myerrmsg("Invalid heap operation, please try again.\n\n");
      }
    }                                 

    switch (type) {
      case BASIC:
        if (nrhs != 6) {
          myprintf("\nUsage: nSwaps = heap('Insert', myHeap,\n");
          myprintf("                     biasedNodeWeight,\n");
          myprintf("                     parentResidualTarget,\n");
          myprintf("                     accumulatedSymbolDecisions,\n");
          myprintf("                     residualDimension);\n\n");
          myprintf("       Inserts a node into the heap having the specified properties.\n");      
          myprintf("       Returns the number of swaps/bubbleUp operations performed.\n\n");
          myerrmsg("Invalid heap operation, please try again.\n\n");
        }                  

        ns = heapInsertB(key, pY, pZ, dim, pH);
        break;

      case LATTICE:
        if (nrhs != 10) {
          myprintf("\nUsage: nSwaps = heap('Insert', myHeap,\n");
          myprintf("                     biasedNodeWeight,\n");
          myprintf("                     parentResidualTarget,\n");
          myprintf("                     accumulatedSymbolDecisions,\n");
          myprintf("                     residualDimension,\n");
          myprintf("                     childOrdinal,\n");
          myprintf("                     unbiasedParentNodeWeight,\n");
          myprintf("                     unconstrainedTarget,\n");
          myprintf("                     previousSiblingDecision);\n\n");
          myprintf("       Inserts a node into the heap having the specified properties.\n");      
          myprintf("       Returns the number of swaps/bubbleUp operations performed.\n\n");
          myerrmsg("Invalid heap operation, please try again.\n\n");
        }                  

        cnum    = (int)mxGetScalar(prhs[6]);
        pweight = (double)mxGetScalar(prhs[7]);
        utarget = (double)mxGetScalar(prhs[8]);
        delta   = (double)mxGetScalar(prhs[9]);
      
        ns = heapInsertI(key, pY, pZ, dim, cnum, pweight, utarget, delta, pH);
        break;

      case LAST:
        if (nrhs != 15) {
          myprintf("\nUsage: nSwaps = heap('Insert', myHeap,\n");
          myprintf("                     biasedNodeWeight,\n");
          myprintf("                     parentResidualTarget,\n");
          myprintf("                     accumulatedSymbolDecisions,\n");
          myprintf("                     residualDimension,\n");
          myprintf("                     childOrdinal,\n");
          myprintf("                     unbiasedParentNodeWeight,\n");
          myprintf("                     unconstrainedTarget,\n");
          myprintf("                     previousSiblingDecision,\n");
          myprintf("                     parentShapingSphereOffset,\n");
          myprintf("                     shapingSphereSquaredRadius,\n");
          myprintf("                     siblingRangeLowerBound,\n");
          myprintf("                     siblingRangeUpperBound,\n");
          myprintf("                     unconstrainedOffset);\n\n");
          myprintf("       Inserts a node into the heap having the specified properties.\n");      
          myprintf("       Returns the number of swaps/bubbleUp operations performed.\n\n");
          myerrmsg("Invalid heap operation, please try again.\n\n");
        } 

        cnum     = (int)mxGetScalar(prhs[6]);
        pweight  = (double)mxGetScalar(prhs[7]);
        utarget  = (double)mxGetScalar(prhs[8]);
        delta    = (double)mxGetScalar(prhs[9]);
        pA       = mxGetPr(prhs[10]);
        sqradius = (double)mxGetScalar(prhs[11]);
        lbound   = (int)mxGetScalar(prhs[12]);