sics/hkl.c

/*---------------------------------------------------------------------------
                         H K L
                         
   Implementation of the SICS object for doing four circle diffractometer
   angle calculations. The actual transformation routines were recoded
   in C from F77 routines provided by J. Allibon, ILL with the MAD
   system.
   
   copyright: see copyright.h
   
   Mark Koennecke, February 1998
   
   Updated to use fourlib.

   Mark Koennecke, December 2001

   Introduced HM mode in order to cope with with the fact that TRICS has 
   three detectors.

   Mark Koennecke, May 2002

   Added handling of the chi ==0 or chi == 180 degree case to tryTweakOmega

   Mark Koennecke, December 2003
   
   Heavily reworked to fit into the new four circle setup
   
   Mark Koennecke, July 2008

   Added calculation of angles between reflection

   Mark Koennecke, March 2013
-----------------------------------------------------------------------------*/
#include <math.h>
#include <ctype.h>
#include "sics.h"
#include "fortify.h"
#include "motor.h"
#include "selector.h"
#include "selvar.h"
#include "fourlib.h"
#include "matrix/matrix.h"
#include "hkl.h"
#include "hkl.i"
#include "splitter.h"
#include "singlex.h"
#include "motorlist.h"
#include "cell.h"
#include "ubfour.h"
/*
  the tolerance in chi we give before we allow to fix omega with phi
*/
#define CHITOLERANCE 3.
#define ABS(x) (x < 0 ? -(x) : (x))
/*-------------------------------------------------------------------------*/
static int HKLSave(void *pData, char *name, FILE * fd)
{
  pHKL self = NULL;

  self = (pHKL) pData;
  if ((self == NULL) || (fd == NULL)) {
    return 1;
  }
  fprintf(fd, "#HKL Settings\n");
  fprintf(fd, "%s scantolerance %f\n", name, self->scanTolerance);
  return 1;
}
/*---------------------------------------------------------------------------*/
static void *HKLInterface(void *pData, int iID)
{
	pHKL self = (pHKL)pData;
	assert(self != NULL);
	
	if(iID == DRIVEID){
		return self->pMotDriv;
	}
	return NULL;
}
/*---------------------------------------------------------------------------*/
static int HKLHalt(void *pData)
{
	pHKL self = (pHKL)pData;
	pSingleDiff single = NULL;
	assert(self != NULL);

	single = SXGetDiffractometer();
	assert(single != NULL);

	return self->pMotList->Halt(&single->motList);
}
/*---------------------------------------------------------------------------*/
static int HKLLimits(void *pData, float val, char *error, int errLen)
{
	pHKL self = (pHKL)pData;
	pSingleDiff single = NULL;
	assert(self != NULL);

	single = SXGetDiffractometer();
	assert(single != NULL);

	return self->pMotList->CheckLimits(&single->motList, val, error, errLen);
}
/*---------------------------------------------------------------------------*/
static long HKLSet(void *pData, SConnection *pCon, float val)
{
	pHKL self = (pHKL)pData;
	pSingleDiff single = NULL;
	assert(self != NULL);

	single = SXGetDiffractometer();
	assert(single != NULL);

	return self->pMotList->SetValue(&single->motList, pCon,val);
}
/*---------------------------------------------------------------------------*/
static int HKLCheck(void *pData, SConnection *pCon)
{
	pHKL self = (pHKL)pData;
	pSingleDiff single = NULL;
	assert(self != NULL);

	single = SXGetDiffractometer();
	assert(single != NULL);

	return self->pMotList->CheckStatus(&single->motList, pCon);
}
/*---------------------------------------------------------------------------*/
static float HKLGet(void *pData, SConnection *pCon)
{
	pHKL self = (pHKL)pData;
	pSingleDiff single = NULL;
	assert(self != NULL);

	single = SXGetDiffractometer();
	assert(single != NULL);

	return self->pMotList->GetValue(&single->motList, pCon);
}
/*---------------------------------------------------------------------------*/
pHKL CreateHKL()
{
  pHKL pNew = NULL;

  /* allocate memory */
  pNew = (pHKL) malloc(sizeof(HKL));
  if (!pNew) {
    return NULL;
  }
  memset(pNew, 0, sizeof(HKL));

  /* create object descriptor */
  pNew->pDes = CreateDescriptor("4-Circle-Calculus");
  pNew->pMotList = makeMotListInterface();
  pNew->pMotDriv = CreateDrivableInterface();
  if (!pNew->pDes || pNew->pMotDriv == NULL || pNew->pMotDriv == NULL) {
    free(pNew);
    return NULL;
  }
  pNew->pDes->SaveStatus = HKLSave;
  pNew->pDes->GetInterface = HKLInterface;
  
  pNew->pMotDriv->Halt = HKLHalt;
  pNew->pMotDriv->CheckLimits = HKLLimits;
  pNew->pMotDriv->CheckStatus = HKLCheck;
  pNew->pMotDriv->SetValue = HKLSet;
  pNew->pMotDriv->GetValue = HKLGet;
  
  pNew->iQuad = 1;
  pNew->iHM = 0;
  pNew->UBinv = NULL;
  pNew->scanTolerance = 2.5;

  return pNew;
}

/*--------------------------------------------------------------------------*/
void DeleteHKL(void *pData)
{
  pHKL self = NULL;

  self = (pHKL) pData;
  if (!self) {
    return;
  }

  if (self->pDes) {
    DeleteDescriptor(self->pDes);
  }
  if (self->pMotDriv) {
    free(self->pMotDriv);
  }

  free(self);
}

/*---------------------------------------------------------------------------*/
int HKLFactory(SConnection * pCon, SicsInterp * pSics, void *pData,
               int argc, char *argv[])
{
  pHKL self = NULL;
  int iRet;

  /* make a new structure */
  self = CreateHKL();
  if (!self) {
    SCWrite(pCon, "ERROR: cannot allocate HKL data structure", eError);
    return 0;
  }

  /* install a command */
  iRet = AddCommand(pSics, "HKL", HKLAction, DeleteHKL, self);
  if (!iRet) {
    SCWrite(pCon, "ERROR: duplicate command HKL not created", eError);
    DeleteHKL(self);
    return 0;
  }
  return 1;
}

/*-------------------------------------------------------------------------*/
void SetHKLScanTolerance(pHKL self, float fVal)
{
  assert(self);
  self->scanTolerance = fVal;
}

/*------------------------------------------------------------------------*/
int GetLambda(pHKL self, float *fVal)
{
  assert(self);
  double val;

  val = SXGetLambda();
  *fVal = (float) val;
  return 1;
}

/*-------------------------------------------------------------------------*/
int GetCurrentHKL(pHKL self, float fHKL[3])
{
  int i;

  assert(self);

  for (i = 0; i < 3; i++) {
    fHKL[i] = self->fLastHKL[i];
  }
  return 1;
}

/*-------------------------------------------------------------------------*/
static void invertMatrix(pHKL self, float fUB[9])
{
  int i;
  MATRIX m;

  if (self->UBinv != NULL) {
    mat_free(self->UBinv);
  }
  m = mat_creat(3, 3, ZERO_MATRIX);
  for (i = 0; i < 3; i++) {
    m[0][i] = fUB[i];
    m[1][i] = fUB[3 + i];
    m[2][i] = fUB[6 + i];
  }
  self->UBinv = mat_inv(m);
  mat_free(m);
}

/*-------------------------------------------------------------------------*/
int SetUB(pHKL self, float fUB[9])
{
  int i;
  double dUB[9];

  assert(self);

  for (i = 0; i < 9; i++) {
    dUB[i] = fUB[i];
  }
  SXSetUB(dUB);
  invertMatrix(self, fUB);

  return 1;
}

/*-------------------------------------------------------------------------*/
int GetUB(pHKL self, float fUB[9])
{
  int i;
  const double *dUB;

  dUB = SXGetUB();
  for (i = 0; i < 9; i++) {
    fUB[i] = (float) dUB[i];
  }
  return 1;
}

/*--------------------------------------------------------------------------*/
int SetNOR(pHKL self, int iNOB)
{
  assert(self);

  /* cannot set normal beam geometry if no nu motor */
  if ((iNOB == 1) && (SXGetMotor(Nu) == NULL))
    return 0;

  switch (iNOB) {
  case 0:
    SXSetMode(Bisecting);
    break;
  case 1:
    SXSetMode(NB);
    break;
  default:
    return 0;
  }
  return 1;
}

  /*-----------------------------------------------------------------------*/
static int checkTheta(pHKL self, double *stt)
{
  char pError[132];
  int iTest;
  float fHard;
  pMotor pTheta;

  pTheta = SXGetMotor(TwoTheta);
  if (pTheta == NULL) {
    return 0;
  }

  iTest = MotorCheckBoundary(pTheta, (float) *stt, &fHard, pError, 131);
  if (!iTest) {
    return -1;
  }
  return iTest;
}

/*----------------------------------------------------------------------*/
int hklInRange(void *data, double fSet[4], int mask[4])
{
  pHKL self = (pHKL) data;
  float fHard, fLimit;
  char pError[132];
  int i, test;
  double dTheta;
  pMotor pOmega, pChi, pPhi;

  pOmega = SXGetMotor(Omega);
  pChi = SXGetMotor(Chi);
  pPhi = SXGetMotor(Phi);
  if (pOmega == NULL || pChi == NULL || pPhi == NULL) {
    return 0;
  }

  /* check two theta  */
  dTheta = fSet[0];
  mask[0] = checkTheta(self, &dTheta);
  fSet[0] = dTheta;

  /* for omega check against the limits +- SCANBORDER in order to allow for
     a omega scan.
   */
  MotorGetPar(pOmega, "softlowerlim", &fLimit);
  if ((float) fSet[1] < fLimit + self->scanTolerance) {
    mask[1] = 1;
  } else {
    mask[1] = 0;
    MotorGetPar(pOmega, "softupperlim", &fLimit);
    if ((float) fSet[1] > fLimit - self->scanTolerance) {
      mask[1] = 0;
    } else {
      mask[1] = 1;
    }
  }

  /* check chi and phi */
  mask[2] = MotorCheckBoundary(pChi, fSet[2], &fHard, pError, 131);
  mask[3] = MotorCheckBoundary(pPhi, fSet[3], &fHard, pError, 131);
  for (i = 0, test = 0; i < 4; i++) {
    test += mask[i];
  }
  if (test != 4) {
    return 0;
  } else {
    return 1;
  }
}

/*-------------------------------------------------------------------------
  calculates the four circle settings. If the position can not be reached
  because of a limit violation, then psi is rotated in 10 degree steps
  until either the loop ends or we finally succeed. If there is a omega
  violation we first try to calculate a delta omega which puts omega
  into the right range. This is a fix because the omega movement is quite
  often restricted due to the cryogenic garbage around the sample.
-------------------------------------------------------------------------*/
int CalculateSettings(pHKL self, float fHKL[3], float fPsi, int iHamil,
                      float fSet[5], SConnection * pCon)
{
  char pBueffel[512];
  double myPsi = fPsi;
  int i, status;
  pSingleDiff single = NULL;
  double dHkl[4], dSet[5];

  /* catch shitty input */
  if ((fHKL[0] == 0.) && (fHKL[1] == 0.) && (fHKL[2] == 0.)) {
    SCWrite(pCon, "ERROR: I will not calculate angles for HKL = (0,0,0) ",
            eError);
    return 0;
  }

  /* some people are stupid.......... */
  myPsi = circlify(fPsi);

  single = SXGetDiffractometer();
  assert(single != NULL);

  for (i = 0; i < 3; i++) {
    dHkl[i] = fHKL[i];
  }
  dHkl[3] = myPsi;

  status = single->calculateSettings(single, dHkl, dSet);
  for (i = 0; i < 5; i++) {
    fSet[i] = dSet[i];
  }

  if (!status) {
    snprintf(pBueffel,sizeof(pBueffel)-1,
            "ERROR: cannot calculate %4.1f %4.1f %4.1f, psi = %4.1f",
            fHKL[0], fHKL[1], fHKL[2], dHkl[3]);
    SCWrite(pCon, pBueffel, eError);
  }
  return status;
}

/*------------------------------------------------------------------------*/
void stopHKLMotors(pHKL self)
{
  pSingleDiff single = NULL;

  single = SXGetDiffractometer();
  assert(single != NULL);
  self->pMotList->Halt(&single->motList);
}

/*------------------------------------------------------------------------*/
int startHKLMotors(pHKL self, SConnection * pCon, float fSet[5])
{
  char pBueffel[512];
  pSingleDiff single = NULL;
  double dSet[5];
  int i, status;

  single = SXGetDiffractometer();
  assert(single != NULL);

  for (i = 0; i < 5; i++) {
    dSet[i] = fSet[i];
  }
  single->settingsToList(single, dSet);

/*  status = self->pMotDriv->SetValue(&single->motList, pCon, .1); */
  status = StartDevice(pServ->pExecutor, "HKL", self->pDes, self, pCon,RUNDRIVE, .1);

  return status;
}

/*-------------------------------------------------------------------------*/
int RunHKL(pHKL self, float fHKL[3],
           float fPsi, int iHamil, SConnection * pCon)
{
  float fSet[5];
  int iRet, i;
  char pBueffel[512];
  pDummy pDum;

  assert(self);
  iRet = CalculateSettings(self, fHKL, fPsi, iHamil, fSet, pCon);
  if (!iRet) {
    SCWrite(pCon, "ERROR: NOT started", eError);
    return 0;
  }

  iRet = startHKLMotors(self, pCon, fSet);
  if (iRet != 1) {
    return iRet;
  }

  for (i = 0; i < 3; i++) {
    self->fLastHKL[i] = fHKL[i];
  }
  return 1;

}

/*-------------------------------------------------------------------------*/
int DriveSettings(pHKL self, float fSet[5], SConnection * pCon)
{
  int iRet;

  iRet = startHKLMotors(self, pCon, fSet);
  if (iRet != 1) {
    stopHKLMotors(self);
  }

  /* wait for end */
  iRet = Wait4Success(pServ->pExecutor);
  switch (iRet) {
  case DEVINT:
    if (SCGetInterrupt(pCon) == eAbortOperation) {
      SCSetInterrupt(pCon, eContinue);
      SCWrite(pCon, "Driving to HKL Aborted", eLogError);
    }
    return 0;
    break;
  case DEVDONE:
    SCWrite(pCon, "Driving to Reflection done", eValue);
    break;
  default:
    SCWrite(pCon, "WARNING: driving to HKL finished with problems",
            eWarning);
    if (SCGetInterrupt(pCon) != eContinue) {
      return 0;
    }
    break;
  }
  return iRet;
}

/*--------------------------------------------------------------------------*/
int DriveHKL(pHKL self, float fHKL[3],
             float fPsi, int iHamil, SConnection * pCon)
{
  int iRet, iReturn;
  long lID;
  char pBueffel[132];

  assert(self);

  /* start running */
  iReturn = 1;
  iRet = RunHKL(self, fHKL, fPsi, iHamil, pCon);
  if (!iRet) {
    StopExe(pServ->pExecutor, "all");
    iReturn = 0;
  }

  /* wait for end */
  iRet = Wait4Success(pServ->pExecutor);
  switch (iRet) {
  case DEVINT:
    if (SCGetInterrupt(pCon) == eAbortOperation) {
      SCSetInterrupt(pCon, eContinue);
      SCWrite(pCon, "Driving to HKL Aborted", eError);
    }
    return 0;
    break;
  case DEVDONE:
    if (fPsi > .01) {
      snprintf(pBueffel, 131,
               "Driving to %8.4f %8.4f %8.4f, psi = %8.4f done",
               fHKL[0], fHKL[1], fHKL[2], fPsi);
    } else {
      snprintf(pBueffel, 131, "Driving to %8.4f %8.4f %8.4f done",
               fHKL[0], fHKL[1], fHKL[2]);
    }
    SCWrite(pCon, pBueffel, eValue);
    break;
  default:
    SCWrite(pCon, "WARNING: driving to HKL finished with problems",
            eWarning);
    if (SCGetInterrupt(pCon) != eContinue) {
      return 0;
    }
    break;
  }
  return iReturn;
}

/*-------------------------------------------------------------------------
  For the conversion from angles to HKL.
  -------------------------------------------------------------------------*/
int GetHKLFromAngles(pHKL self, SConnection * pCon, float fHKL[3])
{
  pSingleDiff single = NULL;
  double dHkl[3];
  int i, status;

  single = SXGetDiffractometer();
  assert(single != NULL);

  status = single->hklFromAngles(single, dHkl);
  for (i = 0; i < 3; i++) {
    fHKL[i] = dHkl[i];
  }
  return status;
}

/*--------------------------------------------------------------------------*/
static int GetCommandData(int argc, char *argv[], float fHKL[3],
                          float *fPsi, int *iHamil, SConnection * pCon)
{
  int iRet, i;
  double d;
  char pBueffel[512];

  if (argc < 3) {
    SCWrite(pCon,
            "ERROR: Insufficient reflection data specified for calculation",
            eError);
    return 0;
  }

  /* get the reflection */
  for (i = 0; i < 3; i++) {
    if (!isNumeric(argv[i])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s is NOT recognized as a number",
              argv[i]);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    fHKL[i] = atof(argv[i]);
  }

  *fPsi = 0.;
  *iHamil = 0;

  /* has psi been specicifed ? */
  if (argc > 3) {
    if (!isNumeric(argv[3])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s is NOT recognized as a number",
              argv[3]);
      SCWrite(pCon, pBueffel, eError);
    }
    *fPsi = atof(argv[3]);
  }

  /* has iHamil been specified ? */
  if (argc > 4) {
    if (!isNumeric(argv[4])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s is NOT recognized as a number",
              argv[4]);
      SCWrite(pCon, pBueffel, eError);
    }
    *iHamil = atof(argv[4]);
  }
  return 1;
}

/*--------------------------------------------------------------------------*/
static int HKLCalculateTheta(pHKL self, float fHKL[3], double *stt)
{
  int status, i;
  MATRIX B, HC, H;
  double theta, d;
  const double *cell;
  reflection r;
  lattice direct;

  B = mat_creat(3, 3, UNIT_MATRIX);
  H = mat_creat(3, 1, ZERO_MATRIX);
  if (B == NULL) {
    return UBNOMEMORY;
  }
  cell = SXGetCell();
  direct.a = cell[0];
  direct.b = cell[1];
  direct.c = cell[2];
  direct.alpha = cell[3];
  direct.beta = cell[4];
  direct.gamma = cell[5];

  status = calculateBMatrix(direct, B);
  if (status < 0) {
    return status;
  }
  for (i = 0; i < 3; i++) {
    H[i][0] = fHKL[i];
  }
  HC = mat_mul(B, H);
  status = calcTheta(SXGetLambda(), HC, &d, &theta);
  mat_free(HC);
  mat_free(B);
  mat_free(H);
  *stt = 2. * theta;
  return status;
}
/*-------------------------------------------------------------------------*/
static int HKLAngle(SConnection *pCon, int argc, char *argv[])
{
  reflection r1, r2;
  lattice direct;
  const double *cell;
  MATRIX B;
  double ang;
  
  if(argc < 8){
    SCWrite(pCon,"ERROR: insufficient no of arguments to hkl angle",eError);
    return 0;
  }
  r1.h = atof(argv[2]);
  r1.k = atof(argv[3]);
  r1.l = atof(argv[4]);

  r2.h = atof(argv[5]);
  r2.k = atof(argv[6]);
  r2.l = atof(argv[7]);
  
  cell = SXGetCell();
  direct.a = cell[0];
  direct.b = cell[1];
  direct.c = cell[2];
  direct.alpha = cell[3];
  direct.beta = cell[4];
  direct.gamma = cell[5];

  B = mat_creat(3,3,ZERO_MATRIX);
  if(B == NULL){
    SCWrite(pCon,"ERROR: out of memory in HKL angle",eError);
    return 0;
  }

  calculateBMatrix(direct,B);
  ang = angleBetweenReflections(B,r1,r2);
  SCPrintf(pCon,eValue,"angle = %f", ang);
  mat_free(B);
  return 1;
}
/*--------------------------------------------------------------------------*/
static int HandleBiToNB(SConnection *pCon, int argc, char *argv[])
{
	double stt, om, chi, phi, gamma, omnb, nu;
	int status;
	
	if(argc < 6){
		SCWrite(pCon,"ERROR: missing arguments, usage: hkl bitonb stt om chi phi", eError);
		return 0;
	}
	status = Tcl_GetDouble(InterpGetTcl(pServ->pSics),argv[2], &stt);
	if(status != 0){
		SCPrintf(pCon,eError,"ERROR: failed to convert %s to double", argv[2]);
		return 0;
	}
	status = Tcl_GetDouble(InterpGetTcl(pServ->pSics),argv[3], &om);
	if(status != 0){
		SCPrintf(pCon,eError,"ERROR: failed to convert %s to double", argv[3]);
		return 0;
	}
	status = Tcl_GetDouble(InterpGetTcl(pServ->pSics),argv[4], &chi);
	if(status != 0){
		SCPrintf(pCon,eError,"ERROR: failed to convert %s to double", argv[4]);
		return 0;
	}
	status = Tcl_GetDouble(InterpGetTcl(pServ->pSics),argv[5], &phi);
	if(status != 0){
		SCPrintf(pCon,eError,"ERROR: failed to convert %s to double", argv[5]);
		return 0;
	}
	status = bisToNormalBeam(stt,om,chi,phi,&omnb, &gamma, &nu);
	if(status != 1){
		SCWrite(pCon,"ERROR: not in range for normal beam", eError);
		return 0;
	}
	SCPrintf(pCon,eValue, "NB = %f,%f,%f", gamma, omnb, nu);
	return 1;
}
/*--------------------------------------------------------------------------*/
int HKLAction(SConnection * pCon, SicsInterp * pSics, void *pData,
              int argc, char *argv[])
{
  int iRet, i, iHamil;
  char pBueffel[512];
  float fUB[9], fPsi, fVal;
  float fHKL[3], fSet[5];
  double dHKL[3], dSet[5];
  double dVal, lambda, stt;
  const double *dUB;
  pHKL self = NULL;
  CommandList *pCom = NULL;
  pDummy pDum = NULL;
  pSingleDiff single = NULL;
  int mode, nargs;

  assert(pCon);
  assert(pSics);

  self = (pHKL) pData;
  assert(self);

  single = SXGetDiffractometer();

  /* enough arguments ? */
  if (argc < 2) {
    SCWrite(pCon, "Insufficient number of arguments to HKL", eError);
    return 0;
  }

/*-------- handle list */
  strtolower(argv[1]);
  if (strcmp(argv[1], "list") == 0) {
    snprintf(pBueffel,sizeof(pBueffel)-1,
            "lambda = %f   Normal Beam = %d  Quadrant = %d HM = %d",
            SXGetLambda(), self->iNOR, self->iQuad, self->iHM);
    SCWrite(pCon, pBueffel, eValue);
    dUB = SXGetUB();
    snprintf(pBueffel,sizeof(pBueffel)-1, "UB = { %f %f %f", dUB[0], dUB[1], dUB[2]);
    SCWrite(pCon, pBueffel, eValue);
    snprintf(pBueffel,sizeof(pBueffel)-1, "       %f %f %f", dUB[3], dUB[4], dUB[5]);
    SCWrite(pCon, pBueffel, eValue);
    snprintf(pBueffel,sizeof(pBueffel)-1, "       %f %f %f }", dUB[6], dUB[7], dUB[8]);
    SCWrite(pCon, pBueffel, eValue);
    return 1;
    snprintf(pBueffel,sizeof(pBueffel)-1, "Last HKL:  %f %f %f ",
            self->fLastHKL[0], self->fLastHKL[1], self->fLastHKL[2]);
    SCWrite(pCon, pBueffel, eValue);
    snprintf(pBueffel, 510, "%s.scantolerance = %f", argv[0],
             self->scanTolerance);
    SCWrite(pCon, pBueffel, eValue);
    return 1;
  }
/*----------- current */
  else if (strcmp(argv[1], "current") == 0) {
    GetHKLFromAngles(self, pCon, fHKL);
    snprintf(pBueffel,sizeof(pBueffel)-1, "Current HKL:  %8.4f %8.4f %8.4f ",
            fHKL[0], fHKL[1], fHKL[2]);
    SCWrite(pCon, pBueffel, eValue);
    return 1;
  } else if (strcmp(argv[1], "fromangles") == 0) {
	mode = SXGetMode();
	if(mode == BiNB) {
		nargs = 7;
	} else  {
		nargs = 6;
	}
    if (argc < nargs) {
      SCWrite(pCon,
              "ERROR: need stt, om, chi,phi to calculate HKL from angles",
              eError);
      return 0;
    }
    for (i = 0; i < nargs-2; i++) {
      iRet = Tcl_GetDouble(InterpGetTcl(pSics), argv[i + 2], &dVal);
      if (iRet != TCL_OK) {
        snprintf(pBueffel, 511, "ERROR: failed to convert %s to number",
                 argv[i + 2]);
        SCWrite(pCon, pBueffel, eError);
        return 0;
      }
      dSet[i] = dVal;
    }
    single->hklFromAnglesGiven(single, dSet, dHKL);
    for (i = 0; i < 3; i++) {
      fHKL[i] = dHKL[i];
    }
    snprintf(pBueffel,sizeof(pBueffel)-1, "HKL from angles:  %8.4f %8.4f %8.4f ",
            fHKL[0], fHKL[1], fHKL[2]);
    SCWrite(pCon, pBueffel, eValue);
    return 1;
  }
/*------------- lambda */
  else if (strcmp(argv[1], "lambda") == 0) {
    lambda = SXGetLambda();
    snprintf(pBueffel, 132, "%s.lambda = %f", argv[0], lambda);
    SCWrite(pCon, pBueffel, eValue);
  }
/*------------- getub*/
  else if (strcmp(argv[1], "getub") == 0) {
    dUB = SXGetUB();
    snprintf(pBueffel, 510, "%s.ub = %f %f %f %f %f %f %f %f %f",
             argv[0], dUB[0], dUB[1], dUB[2],
             dUB[3], dUB[4], dUB[5], dUB[6], dUB[7], dUB[8]);
    SCWrite(pCon, pBueffel, eValue);
    return 1;
  }
/*------------ UB */
  else if (strcmp(argv[1], "setub") == 0) {
    if (argc < 11) {
      SCWrite(pCon, "ERROR: Insufficient number of arguments to HKL setUB",
              eError);
      return 0;
    }
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    for (i = 2; i < 11; i++) {
      if (!isNumeric(argv[i])) {
        snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s was not recognized as a number",
                argv[i]);
        SCWrite(pCon, pBueffel, eError);
        return 0;
      }
      fUB[i - 2] = atof(argv[i]);
    }
    iRet = SetUB(self, fUB);
    if (!iRet) {
      SCWrite(pCon, "ERROR: cannot set UB Matrix", eError);
      return 0;
    }
    SCSendOK(pCon);
    return 1;
  }
/*------------- HM mode */
  else if (strcmp(argv[1], "hm") == 0) {
    if (argc < 3) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "%s.hm = %d", argv[0], self->iHM);
      SCWrite(pCon, pBueffel, eValue);
      return 1;
    }
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    if (!isNumeric(argv[2])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s was not recognized as a number",
              argv[2]);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    self->iHM = atoi(argv[2]);
    SCSendOK(pCon);
    return 1;
  }
/*------------- normal beam */
  else if (strcmp(argv[1], "nb") == 0) {
    if (argc < 3) {
      snprintf(pBueffel, 511, "%s.nb = %d", argv[0], self->iNOR);
      SCWrite(pCon, pBueffel, eValue);
      return 1;
    }
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    if (!isNumeric(argv[2])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s was not recognized as a number",
              argv[2]);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    iRet = SetNOR(self, atoi(argv[2]));
    if (!iRet) {
      SCWrite(pCon,
              "ERROR: cannot set Normal Beam geometry, probably nu motor undefined",
              eError);
      return 0;
    }
    SCSendOK(pCon);
    return 1;
  }
/*------------- quadrant  */
  else if (strcmp(argv[1], "quadrant") == 0) {
    if (argc < 3) {
      SCWrite(pCon,
              "ERROR: Insufficient number of arguments to HKL quadrant",
              eError);
      return 0;
    }
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    if (!isNumeric(argv[2])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s was not recognized as a number",
              argv[2]);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    iRet = atoi(argv[2]);
    if (!((iRet == 1) || (iRet == 0))) {
      snprintf(pBueffel,sizeof(pBueffel)-1,
              "ERROR: Invalid parameter %d for quadrant, posiible: 0,1",
              iRet);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    self->iQuad = iRet;
    SCSendOK(pCon);
    return 1;
  }
/*------------- scantolerance  */
  else if (strcmp(argv[1], "scantolerance") == 0) {
    if (argc < 3) {
      snprintf(pBueffel, 510, "%s.scantolerance = %f", argv[0],
               self->scanTolerance);
      SCWrite(pCon, pBueffel, eValue);
      return 1;
    }
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    if (!isNumeric(argv[2])) {
      snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: %s was not recognized as a number",
              argv[2]);
      SCWrite(pCon, pBueffel, eError);
      return 0;
    }
    self->scanTolerance = atof(argv[2]);
    SCSendOK(pCon);
    return 1;
  }
/*------------- calculate */
  else if (strcmp(argv[1], "calc") == 0) {
    iRet = GetCommandData(argc - 2, &argv[2], fHKL, &fPsi, &iHamil, pCon);
    if (!iRet) {
      return 0;
    }
    iRet = CalculateSettings(self, fHKL, fPsi, iHamil, fSet, pCon);
    if (SXGetMode() == NB) {
      snprintf(pBueffel,sizeof(pBueffel)-1, " gamma = %f, omega = %f, nu = %f",
              fSet[0], fSet[1], fSet[2]);
      SCWrite(pCon, pBueffel, eValue);
    } else if(SXGetMode() == BiNB){
      snprintf(pBueffel,sizeof(pBueffel)-1, " 2-theta = %f, omega = %f, chi = %f, phi = %f, nu = %f",
              fSet[0], fSet[1], fSet[2], fSet[3], fSet[4]);
      SCWrite(pCon, pBueffel, eValue);
    } else {
      snprintf(pBueffel,sizeof(pBueffel)-1, " 2-theta = %f, omega = %f, chi = %f, phi = %f",
              fSet[0], fSet[1], fSet[2], fSet[3]);
      SCWrite(pCon, pBueffel, eValue);
    }
    if (!iRet) {
      SCWrite(pCon,
              "WARNING: Cannot drive to the hkl of your desire", eWarning);
      return 0;
    }
    return 1;
  }
/*------------- calculate theta*/
  else if (strcmp(argv[1], "calctth") == 0) {
    iRet = GetCommandData(argc - 2, &argv[2], fHKL, &fPsi, &iHamil, pCon);
    if (!iRet) {
      return 0;
    }
    HKLCalculateTheta(self, fHKL, &stt);
    SCPrintf(pCon, eValue, "two-theta = %lf", stt);
    return 1;
  }
  /*--------------- calculate angle */
  else if (strcmp(argv[1], "angle") == 0) {
    return HKLAngle(pCon,argc,argv);
  }
/*------------------ run */
  else if (strcmp(argv[1], "run") == 0) {
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    iRet = GetCommandData(argc - 2, &argv[2], fHKL, &fPsi, &iHamil, pCon);
    if (!iRet) {
      return 0;
    }
    iRet = RunHKL(self, fHKL, fPsi, iHamil, pCon);
    if (!iRet) {
      return 0;
    } else {
      SCSendOK(pCon);
      return 1;
    }
  }
/*------------------ drive */
  else if (strcmp(argv[1], "drive") == 0) {
    if (!SCMatchRights(pCon, usUser)) {
      return 0;
    }
    iRet = GetCommandData(argc - 2, &argv[2], fHKL, &fPsi, &iHamil, pCon);
    if (!iRet) {
      return 0;
    }
    iRet = DriveHKL(self, fHKL, fPsi, iHamil, pCon);
    if (!iRet) {
      return 0;
    } else {
      SCSendOK(pCon);
      return 1;
    }
  }
    /*------------------ bitonb */
    else if (strcmp(argv[1], "bitonb") == 0) {
    	return HandleBiToNB(pCon, argc, argv);
  } else {
    snprintf(pBueffel,sizeof(pBueffel)-1, "ERROR: subcommand --> %s <-- not recognized",
            argv[1]);
    SCWrite(pCon, pBueffel, eError);
    return 0;
  }
  return 0;                     /* not reached */
}