1 /*
2   Space Group Info's (c) 1994-96 Ralf W. Grosse-Kunstleve
3  */
4 
5 #include <stdio.h>
6 #include <stdlib.h>
7 
8 
9 #include "sginfo.h"
10 
11 
12 static const char *IErr_Inc_SymMx =
13   "Internal Error: Inconsistent symmetry matrices";
14 
15 
IsSysAbsent_hkl(const T_SgInfo * SgInfo,int h,int k,int l,int * TH_Restriction)16 int IsSysAbsent_hkl(const T_SgInfo *SgInfo,
17                     int h, int k, int l, int *TH_Restriction)
18 {
19   int           iTrV, nTrV;
20   const int     *TrV;
21   int           iList, mh, mk, ml, hm, km, lm;
22   int           TH, THr, FlagMismatch;
23   const T_RTMx  *lsmx;
24 
25 
26   mh = -h;
27   mk = -k;
28   ml = -l;
29 
30   /* check list of symmetry operations
31      take care of lattice type and "centric" flag */
32 
33   THr = -1;
34   if (TH_Restriction != NULL) *TH_Restriction = THr;
35   FlagMismatch = 0;
36 
37   nTrV = SgInfo->LatticeInfo->nTrVector;
38   lsmx = SgInfo->ListSeitzMx;
39 
40   for (iList = 0; iList < SgInfo->nList; iList++, lsmx++)
41   {
42     hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
43     km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
44     lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
45 
46     TrV = SgInfo->LatticeInfo->TrVector;
47 
48     for (iTrV = 0; iTrV < nTrV; iTrV++)
49     {
50       TH =  (lsmx->s.T[0] + *TrV++) * h;
51       TH += (lsmx->s.T[1] + *TrV++) * k;
52       TH += (lsmx->s.T[2] + *TrV++) * l;
53       TH %= STBF; if (TH < 0) TH += STBF;
54 
55       if      (mh == hm && mk == km && ml == lm)
56       {
57         if (TH != 0 && SgInfo->Centric == -1)
58           return -(iList + 1 + iTrV * SgInfo->nList);
59 
60         if (THr < 0)
61           THr = TH;
62         else if (THr != TH)
63           FlagMismatch = 1; /* must be systematic absent */
64                             /* will check later ...      */
65       }
66       else if ( h == hm &&  k == km &&  l == lm)
67       {
68         if (TH != 0)
69           return  (iList + 1 + iTrV * SgInfo->nList);
70       }
71       else
72         break;
73     }
74   }
75 
76   if (THr >= 0 && FlagMismatch) /* ... consistency check */
77     SetSgError(IErr_Inc_SymMx);
78 
79   if (TH_Restriction != NULL)
80   {
81     if (SgInfo->Centric == -1) *TH_Restriction = 0;
82     else                       *TH_Restriction = THr;
83   }
84 
85   return 0;
86 }
87 
88 
BuildEq_hkl(const T_SgInfo * SgInfo,T_Eq_hkl * Eq_hkl,int h,int k,int l)89 int BuildEq_hkl(const T_SgInfo *SgInfo, T_Eq_hkl *Eq_hkl, int h, int k, int l)
90 {
91   int       iList, hm, km, lm, i;
92   T_RTMx    *lsmx;
93   T_Eq_hkl  BufEq_hkl;
94 
95 
96   if (Eq_hkl == NULL)
97     Eq_hkl = &BufEq_hkl;
98 
99   Eq_hkl->M = 1;
100   Eq_hkl->N = 1;
101   Eq_hkl->h[0] = h;
102   Eq_hkl->k[0] = k;
103   Eq_hkl->l[0] = l;
104   Eq_hkl->TH[0] = 0;
105 
106   if (! (h || k || l))
107     return Eq_hkl->M; /* this is 000 */
108 
109   Eq_hkl->M++;
110 
111   /* check list of symmetry operations */
112 
113   lsmx = &SgInfo->ListSeitzMx[1]; /* skip first = identity matrix */
114 
115   for (iList = 1; iList < SgInfo->nList; iList++, lsmx++)
116   {
117     hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
118     km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
119     lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
120 
121     for (i = 0; i < Eq_hkl->N; i++)
122     {
123       if ( ( hm == Eq_hkl->h[i] &&  km == Eq_hkl->k[i] &&  lm == Eq_hkl->l[i])
124         || (-hm == Eq_hkl->h[i] && -km == Eq_hkl->k[i] && -lm == Eq_hkl->l[i]))
125         break;
126     }
127 
128     if (i == Eq_hkl->N)
129     {
130       if (Eq_hkl->N >= 24) {
131         SetSgError(IErr_Inc_SymMx);
132         return 0;
133       }
134 
135       Eq_hkl->h[i] = hm;
136       Eq_hkl->k[i] = km;
137       Eq_hkl->l[i] = lm;
138 
139           Eq_hkl->TH[i] = (  lsmx->s.T[0] * h
140                            + lsmx->s.T[1] * k
141                            + lsmx->s.T[2] * l) % STBF;
142       if (Eq_hkl->TH[i] < 0)
143           Eq_hkl->TH[i] += STBF;
144 
145       Eq_hkl->M += 2;
146       Eq_hkl->N++;
147     }
148   }
149 
150   if (SgInfo->nList % Eq_hkl->N) /* another error trap */ {
151     SetSgError(IErr_Inc_SymMx);
152     return 0;
153   }
154 
155   return Eq_hkl->M;
156 }
157 
158 
AreSymEquivalent_hkl(const T_SgInfo * SgInfo,int h1,int k1,int l1,int h2,int k2,int l2)159 int AreSymEquivalent_hkl(const T_SgInfo *SgInfo, int h1, int k1, int l1,
160                                                  int h2, int k2, int l2)
161 {
162   int     iList, mh2, mk2, ml2, hm, km, lm;
163   T_RTMx  *lsmx;
164 
165 
166   mh2 = -h2;
167   mk2 = -k2;
168   ml2 = -l2;
169 
170   /* check list of symmetry operations */
171 
172   lsmx = SgInfo->ListSeitzMx;
173 
174   for (iList = 0; iList < SgInfo->nList; iList++, lsmx++)
175   {
176     hm = lsmx->s.R[0] * h1 + lsmx->s.R[3] * k1 + lsmx->s.R[6] * l1;
177     km = lsmx->s.R[1] * h1 + lsmx->s.R[4] * k1 + lsmx->s.R[7] * l1;
178     lm = lsmx->s.R[2] * h1 + lsmx->s.R[5] * k1 + lsmx->s.R[8] * l1;
179 
180     if      ( h2 == hm &&  k2 == km &&  l2 == lm)
181       return  (iList + 1);
182 
183     else if (mh2 == hm && mk2 == km && ml2 == lm)
184       return -(iList + 1);
185   }
186 
187   return 0;
188 }
189 
190 
SetListMin_hkl(const T_SgInfo * SgInfo,int Maxk,int Maxl,int * Minh,int * Mink,int * Minl)191 void SetListMin_hkl(const T_SgInfo *SgInfo,            int  Maxk, int  Maxl,
192                                             int *Minh, int *Mink, int *Minl)
193 {
194   *Minh = 0;
195 
196   switch(SgInfo->XtalSystem)
197   {
198     case XS_Triclinic:
199       *Mink = -Maxk;
200       *Minl = -Maxl;
201       break;
202     case XS_Monoclinic:
203       if (SgInfo->UniqueRefAxis == 'z')
204       {
205         *Mink = -Maxk;
206         *Minl = 0;
207       }
208       else
209       {
210         *Mink = 0;
211         *Minl = -Maxl;
212       }
213       break;
214     default:
215       if (SgInfo->XtalSystem == XS_Trigonal && SgInfo->UniqueDirCode == '*')
216         *Mink = -Maxk;
217       else
218         *Mink = 0;
219       *Minl = 0;
220       break;
221   }
222 }
223 
224 
IsSuppressed_hkl(const T_SgInfo * SgInfo,int Minh,int Mink,int Minl,int Maxk,int Maxl,int h,int k,int l)225 int IsSuppressed_hkl(const T_SgInfo *SgInfo, int Minh, int Mink, int Minl,
226                                                        int Maxk, int Maxl,
227                                              int    h, int    k, int    l)
228 {
229   int     iList, mate, hm, km, lm;
230   T_RTMx  *lsmx;
231 
232 
233   /* check for Friedel mate first */
234 
235   hm = -h, km = -k, lm = -l;
236 
237   if (   (Minh <= hm && hm <=    h)
238       && (Mink <= km && km <= Maxk)
239       && (Minl <= lm && lm <= Maxl))
240   {
241     if (hm < h) return -1;
242     else /* if (h == 0) */
243       if (km < k) return -1;
244       else if (k == 0)
245         if (lm < l) return -1;
246   }
247   lsmx = &SgInfo->ListSeitzMx[1]; /* skip first = identity matrix */
248 
249   for (iList = 1; iList < SgInfo->nList; iList++, lsmx++)
250   {
251     /* check if equivalent hm, km, lm are inside loop range ... */
252 
253     hm = lsmx->s.R[0] * h + lsmx->s.R[3] * k + lsmx->s.R[6] * l;
254     km = lsmx->s.R[1] * h + lsmx->s.R[4] * k + lsmx->s.R[7] * l;
255     lm = lsmx->s.R[2] * h + lsmx->s.R[5] * k + lsmx->s.R[8] * l;
256 
257     for (mate = 0; mate < 2; mate++)
258     {
259       if (mate) hm = -hm, km = -km, lm = -lm; /* ... or friedel mate */
260 
261       if (   Minh <= hm && hm <=    h
262           && Mink <= km && km <= Maxk
263           && Minl <= lm && lm <= Maxl)
264       {
265         /* ... and were processed before */
266 
267         if (hm < h)
268           return (mate ? -(iList + 1) : iList + 1);
269         else /* if (hm == h) */
270           if (km < k)
271             return (mate ? -(iList + 1) : iList + 1);
272           else if (km == k)
273             if (lm < l)
274               return (mate ? -(iList + 1) : iList + 1);
275       }
276     }
277   }
278 
279   return 0;
280 }
281