Files
sics/sghkl.c
koennecke 91d4af0541 - Adapted indenation to new agreed upon system
- Added support for second generation scriptcontext based counter
2009-02-13 09:00:03 +00:00

269 lines
6.2 KiB
C

/*
Space Group Info's (c) 1994-96 Ralf W. Grosse-Kunstleve
*/
#include <stdio.h>
#include <stdlib.h>
#include "sginfo.h"
static const char *IErr_Inc_SymMx =
"Internal Error: Inconsistent symmetry matrices";
int IsSysAbsent_hkl(const T_SgInfo * SgInfo,
int h, int k, int l, int *TH_Restriction)
{
int iTrV, nTrV;
const int *TrV;
int iList, mh, mk, ml, hm, km, lm;
int TH, THr, FlagMismatch;
const T_RTMx *lsmx;
mh = -h;
mk = -k;
ml = -l;
/* check list of symmetry operations
take care of lattice type and "centric" flag */
THr = -1;
if (TH_Restriction != NULL)
*TH_Restriction = THr;
FlagMismatch = 0;
nTrV = SgInfo->LatticeInfo->nTrVector;
lsmx = SgInfo->ListSeitzMx;
for (iList = 0; iList < SgInfo->nList; iList++, lsmx++) {
hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
TrV = SgInfo->LatticeInfo->TrVector;
for (iTrV = 0; iTrV < nTrV; iTrV++) {
TH = (lsmx->s.T[0] + *TrV++) * h;
TH += (lsmx->s.T[1] + *TrV++) * k;
TH += (lsmx->s.T[2] + *TrV++) * l;
TH %= STBF;
if (TH < 0)
TH += STBF;
if (mh == hm && mk == km && ml == lm) {
if (TH != 0 && SgInfo->Centric == -1)
return -(iList + 1 + iTrV * SgInfo->nList);
if (THr < 0)
THr = TH;
else if (THr != TH)
FlagMismatch = 1; /* must be systematic absent */
/* will check later ... */
} else if (h == hm && k == km && l == lm) {
if (TH != 0)
return (iList + 1 + iTrV * SgInfo->nList);
} else
break;
}
}
if (THr >= 0 && FlagMismatch) /* ... consistency check */
SetSgError(IErr_Inc_SymMx);
if (TH_Restriction != NULL) {
if (SgInfo->Centric == -1)
*TH_Restriction = 0;
else
*TH_Restriction = THr;
}
return 0;
}
int BuildEq_hkl(const T_SgInfo * SgInfo, T_Eq_hkl * Eq_hkl, int h, int k,
int l)
{
int iList, hm, km, lm, i;
T_RTMx *lsmx;
T_Eq_hkl BufEq_hkl;
if (Eq_hkl == NULL)
Eq_hkl = &BufEq_hkl;
Eq_hkl->M = 1;
Eq_hkl->N = 1;
Eq_hkl->h[0] = h;
Eq_hkl->k[0] = k;
Eq_hkl->l[0] = l;
Eq_hkl->TH[0] = 0;
if (!(h || k || l))
return Eq_hkl->M; /* this is 000 */
Eq_hkl->M++;
/* check list of symmetry operations */
lsmx = &SgInfo->ListSeitzMx[1]; /* skip first = identity matrix */
for (iList = 1; iList < SgInfo->nList; iList++, lsmx++) {
hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
for (i = 0; i < Eq_hkl->N; i++) {
if ((hm == Eq_hkl->h[i] && km == Eq_hkl->k[i] && lm == Eq_hkl->l[i])
|| (-hm == Eq_hkl->h[i] && -km == Eq_hkl->k[i]
&& -lm == Eq_hkl->l[i]))
break;
}
if (i == Eq_hkl->N) {
if (Eq_hkl->N >= 24) {
SetSgError(IErr_Inc_SymMx);
return 0;
}
Eq_hkl->h[i] = hm;
Eq_hkl->k[i] = km;
Eq_hkl->l[i] = lm;
Eq_hkl->TH[i] = (lsmx->s.T[0] * h
+ lsmx->s.T[1] * k + lsmx->s.T[2] * l) % STBF;
if (Eq_hkl->TH[i] < 0)
Eq_hkl->TH[i] += STBF;
Eq_hkl->M += 2;
Eq_hkl->N++;
}
}
if (SgInfo->nList % Eq_hkl->N) { /* another error trap */
SetSgError(IErr_Inc_SymMx);
return 0;
}
return Eq_hkl->M;
}
int AreSymEquivalent_hkl(const T_SgInfo * SgInfo, int h1, int k1, int l1,
int h2, int k2, int l2)
{
int iList, mh2, mk2, ml2, hm, km, lm;
T_RTMx *lsmx;
mh2 = -h2;
mk2 = -k2;
ml2 = -l2;
/* check list of symmetry operations */
lsmx = SgInfo->ListSeitzMx;
for (iList = 0; iList < SgInfo->nList; iList++, lsmx++) {
hm = lsmx->s.R[0] * h1 + lsmx->s.R[3] * k1 + lsmx->s.R[6] * l1;
km = lsmx->s.R[1] * h1 + lsmx->s.R[4] * k1 + lsmx->s.R[7] * l1;
lm = lsmx->s.R[2] * h1 + lsmx->s.R[5] * k1 + lsmx->s.R[8] * l1;
if (h2 == hm && k2 == km && l2 == lm)
return (iList + 1);
else if (mh2 == hm && mk2 == km && ml2 == lm)
return -(iList + 1);
}
return 0;
}
void SetListMin_hkl(const T_SgInfo * SgInfo, int Maxk, int Maxl,
int *Minh, int *Mink, int *Minl)
{
*Minh = 0;
switch (SgInfo->XtalSystem) {
case XS_Triclinic:
*Mink = -Maxk;
*Minl = -Maxl;
break;
case XS_Monoclinic:
if (SgInfo->UniqueRefAxis == 'z') {
*Mink = -Maxk;
*Minl = 0;
} else {
*Mink = 0;
*Minl = -Maxl;
}
break;
default:
if (SgInfo->XtalSystem == XS_Trigonal && SgInfo->UniqueDirCode == '*')
*Mink = -Maxk;
else
*Mink = 0;
*Minl = 0;
break;
}
}
int IsSuppressed_hkl(const T_SgInfo * SgInfo, int Minh, int Mink, int Minl,
int Maxk, int Maxl, int h, int k, int l)
{
int iList, mate, hm, km, lm;
T_RTMx *lsmx;
/* check for Friedel mate first */
hm = -h, km = -k, lm = -l;
if ((Minh <= hm && hm <= h)
&& (Mink <= km && km <= Maxk)
&& (Minl <= lm && lm <= Maxl)) {
if (hm < h)
return -1;
else /* if (h == 0) */ if (km < k)
return -1;
else if (k == 0)
if (lm < l)
return -1;
}
lsmx = &SgInfo->ListSeitzMx[1]; /* skip first = identity matrix */
for (iList = 1; iList < SgInfo->nList; iList++, lsmx++) {
/* check if equivalent hm, km, lm are inside loop range ... */
hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
for (mate = 0; mate < 2; mate++) {
if (mate)
hm = -hm, km = -km, lm = -lm; /* ... or friedel mate */
if (Minh <= hm && hm <= h
&& Mink <= km && km <= Maxk && Minl <= lm && lm <= Maxl) {
/* ... and were processed before */
if (hm < h)
return (mate ? -(iList + 1) : iList + 1);
else /* if (hm == h) */ if (km < k)
return (mate ? -(iList + 1) : iList + 1);
else if (km == k)
if (lm < l)
return (mate ? -(iList + 1) : iList + 1);
}
}
}
return 0;
}