Develop and Download Open Source Software

Browse Subversion Repository

Contents of /Conograph/trunk/src/SortingLattice.cc

Parent Directory Parent Directory | Revision Log Revision Log


Revision 25 - (show annotations) (download) (as text)
Mon Jul 7 02:35:51 2014 UTC (9 years, 8 months ago) by rtomiyasu
File MIME type: text/x-c++src
File size: 24185 byte(s)
Source codes of version 0.9.99
1 /*
2 * The MIT License
3
4 Conograph (powder auto-indexing program)
5
6 Copyright (c) <2012> <Ryoko Oishi-Tomiyasu, KEK>
7
8 Permission is hereby granted, free of charge, to any person obtaining a copy
9 of this software and associated documentation files (the "Software"), to deal
10 in the Software without restriction, including without limitation the rights
11 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
12 copies of the Software, and to permit persons to whom the Software is
13 furnished to do so, subject to the following conditions:
14
15 The above copyright notice and this permission notice shall be included in
16 all copies or substantial portions of the Software.
17
18 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
19 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
20 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
21 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
22 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
23 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
24 THE SOFTWARE.
25 *
26 */
27 #ifdef _OPENMP
28 # include <omp.h>
29 #endif
30 #include "utility_func/chToDouble.hh"
31 #include "utility_func/transform_sym_matrix.hh"
32 #include "utility_lattice_reduction/put_Selling_reduced_lattice.hh"
33 #include "lattice_symmetry/VCLatticeFigureOfMeritToCheckSymmetry.hh"
34 #include "lattice_symmetry/ReducedLatticeToCheckEquiv.hh"
35 #include "zerror_type/error_out.hh"
36 #include "zlog/zlog.hh"
37 #include "ControlParam.hh"
38 #include "utility_func/stopx.hh"
39 #include "SortingLattice.hh"
40
41 const bool SortingLattice::m_DoesPrudentSymSearch = false;
42 const Double SortingLattice::m_cv2 = 0.5;
43
44 SortingLattice::SortingLattice()
45 {
46 for(Int4 i=0; i<NUM_LS; i++)
47 {
48 OutputSymmetry[i] = false;
49 JudgeSymmetry[i] = false;
50 }
51
52 m_resol = 0.0;
53 m_num_ref_figure_of_merit = 20;
54 m_etype_peak_shift = kPeakShiftFunction_Type0;
55 m_WlengthX = 1.54056;
56 }
57
58
59 SortingLattice::~SortingLattice()
60 {
61 }
62
63
64 // Set the member variables.
65 void SortingLattice::setParam(const ControlParam& cont)
66 {
67 OutputSymmetry[(size_t)Triclinic] = cont.putOutputSymmetry(Triclinic);
68 JudgeSymmetry[(size_t)Triclinic] = false;
69 for(Int4 i=1; i<NUM_LS; i++)
70 {
71 OutputSymmetry[i] = cont.putOutputSymmetry(eBravaisType(i));
72 JudgeSymmetry[i] = cont.putOutputSymmetry(eBravaisType(i));
73 }
74
75 if( JudgeSymmetry[(size_t)Cubic_P] )
76 {
77 JudgeSymmetry[(size_t)Tetragonal_P] = true;
78 }
79 if( JudgeSymmetry[(size_t)Hexagonal] )
80 {
81 JudgeSymmetry[(size_t)Monoclinic_P] = true;
82 }
83 if( JudgeSymmetry[(size_t)Tetragonal_P] )
84 {
85 JudgeSymmetry[(size_t)Orthorhombic_P] = true;
86 }
87 if( JudgeSymmetry[(size_t)Orthorhombic_P] )
88 {
89 JudgeSymmetry[(size_t)Monoclinic_P] = true;
90 }
91
92 if( JudgeSymmetry[(size_t)Orthorhombic_C] )
93 {
94 JudgeSymmetry[(size_t)Monoclinic_B] = true;
95 }
96
97 if( JudgeSymmetry[(size_t)Cubic_I] )
98 {
99 JudgeSymmetry[(size_t)Tetragonal_I] = true;
100 }
101 if( JudgeSymmetry[(size_t)Tetragonal_I] )
102 {
103 JudgeSymmetry[(size_t)Orthorhombic_I] = true;
104 }
105
106 if( JudgeSymmetry[(size_t)Cubic_F] )
107 {
108 JudgeSymmetry[(size_t)Orthorhombic_F] = true;
109 }
110
111 m_resol = cont.putResolution();
112 m_num_ref_figure_of_merit = cont.putNumberOfReflectionsForFigureOfMerit();
113 m_etype_peak_shift = cont.putPeakShiftFunctionType();
114 m_WlengthX = cont.putWaveLength();
115
116 const vector<Double>& peak_shift_param_rad = cont.putPeakShiftParamRadian();
117 const Int4 param_num = peak_shift_param_rad.size();
118 assert( m_etype_peak_shift != kPeakShiftFunction_Type0 || param_num == 0 );
119 assert( m_etype_peak_shift != kPeakShiftFunction_Type1 || param_num == 1 );
120
121 m_peak_shift_param_rad.resize(param_num);
122 for(Int4 i=0; i<param_num; i++)
123 {
124 m_peak_shift_param_rad[i] = peak_shift_param_rad[i];
125 }
126 }
127
128
129
130
131 void SortingLattice::putCentringTypes(const ReducedVCLatticeToCheckBravais& RLCB,
132 const VCLatticeFigureOfMeritToCheckSymmetry& lattice_original,
133 const BravaisType& brat,
134 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result) const
135 {
136 lattice_result.clear();
137
138 const map< SymMat<VCData>, NRMat<Int4> >& S_red_tray = RLCB.checkCentringType(brat);
139 if( S_red_tray.empty() ) return;
140
141 // The lattice of RLCB has at least the symmetry given by eblat.
142 SymMat<VCData> S_super(4);
143 NRMat<Int4> trans_mat(4,3);
144
145 for(map< SymMat<VCData>, NRMat<Int4> >::const_iterator it=S_red_tray.begin(); it!=S_red_tray.end(); it++)
146 {
147 S_super = transform_sym_matrix(it->second, it->first);
148 trans_mat = identity_matrix<Int4>(4);
149
150 // S_super = trans_mat * it->second * it->first * Transpose(trans_mat * it->second) is Delone reduced.
151 if( !put_Selling_reduced_dim_less_than_4(S_super, trans_mat) )
152 {
153 assert( false );
154 }
155 moveSmallerDiagonalLeftUpper(S_super, trans_mat);
156
157 lattice_result.push_back( VCLatticeFigureOfMeritToCheckSymmetry( brat, SymMat43_VCData(it->first, mprod(trans_mat, it->second) ),
158 lattice_original.putLatticeFigureOfMerit().putPeakShiftFunctionType(),
159 lattice_original.putLatticeFigureOfMerit().putWaveLength(),
160 lattice_original.putLatticeFigureOfMerit().putPeakShiftParamRadian() ) );
161 }
162 }
163
164
165
166
167
168
169 void SortingLattice::putLatticeCandidatesForTriclinic(const vector<SymMat43_VCData>& S_super,
170 const Double& MIN_NormM,
171 const Double& MIN_RevM,
172 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_tri) const
173 {
174 const Int4 num_topo = S_super.size();
175 lattice_result_tri.clear();
176
177 /* 2011.10.19 VIC Tamura */
178 Int4 LOOP_COUNTER = 0;
179
180 #ifdef _OPENMP
181 #pragma omp parallel
182 #endif
183 {
184 vector< VecDat3<Int4> > closest_hkl_tray;
185 vector<bool> is_cal_Q_observed_tray;
186 vector<VCLatticeFigureOfMeritToCheckSymmetry> latFOM_tray;
187
188 #ifdef _OPENMP
189 #pragma omp for
190 #endif
191 for(Int4 n=0; n<num_topo; n++)
192 {
193 /* 2011.10.19 VIC Tamura */
194 SET_PROGRESS(100*(LOOP_COUNTER++)/num_topo, 65, 1); // critical, but works
195 if(IS_CANSELED()) continue;
196
197 VCLatticeFigureOfMeritToCheckSymmetry latFOM(BravaisType( pair<eCentringType, ePointGroup>(Prim, Ci) ), S_super[n],
198 m_etype_peak_shift, m_WlengthX, m_peak_shift_param_rad);
199
200 latFOM.setFigureOfMerit(m_num_ref_figure_of_merit,
201 VCData::putPeakQData(),
202 closest_hkl_tray, is_cal_Q_observed_tray);
203
204 LatticeFigureOfMeritZeroShift latFOM2 = latFOM.putLatticeFigureOfMerit();
205 pair<bool, ZErrorMessage> ans = latFOM2.fitLatticeParameterLinear(VCData::putPeakQData(),
206 closest_hkl_tray, is_cal_Q_observed_tray, false);
207
208 if( ans.first )
209 {
210 assert( latFOM.putLatticeFigureOfMerit().putFiguresOfMerit().putNumberOfReflectionsForFigureOfMerit() > 0 );
211 latFOM2.setFigureOfMerit(latFOM.putLatticeFigureOfMerit().putFiguresOfMerit().putNumberOfReflectionsForFigureOfMerit(),
212 VCData::putPeakQData(),
213 closest_hkl_tray, is_cal_Q_observed_tray);
214 if( LatticeFigureOfMerit::cmpFOMWu( latFOM2, latFOM.putLatticeFigureOfMerit() ) )
215 {
216 latFOM.setLatticeFigureOfMerit(latFOM2);
217 }
218 }
219 const LatticeFigureOfMerit::SetOfFigureOfMerit& setFOM = latFOM.putLatticeFigureOfMerit().putFiguresOfMerit();
220 if( setFOM.putFigureOfMeritWu() < MIN_NormM ) continue;
221 if( setFOM.putReversedFigureOfMerit() < MIN_RevM ) continue;
222
223 #ifdef _OPENMP
224 #pragma omp critical
225 #endif
226 {
227 lattice_result_tri.push_back( latFOM );
228 }
229 }
230 }
231
232 /* 2011.10.19 VIC Tamura */
233 CHECK_INTERRUPTION();
234
235 // sort( lattice_result_tri.begin(), lattice_result_tri.end() );
236 }
237
238
239 void SortingLattice::putLatticeCandidatesForEachBravaisTypes(
240 const Double& MIN_NormM,
241 const Double& MIN_RevM,
242 const Int4& MAX_SIZE,
243 const eABCaxis& abc_axis,
244 const eRHaxis& rh_axis,
245 vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[NUM_LS]) const
246 {
247 try{
248
249 for(Int4 i=1; i<NUM_LS; i++)
250 {
251 lattice_result[i].clear();
252 }
253 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_tri = lattice_result[(size_t)Triclinic];
254 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_mono_P = lattice_result[(size_t)Monoclinic_P];
255 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_mono_B = lattice_result[(size_t)Monoclinic_B];
256 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_ortho_P = lattice_result[(size_t)Orthorhombic_P];
257 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_ortho_B = lattice_result[(size_t)Orthorhombic_C];
258 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_ortho_I = lattice_result[(size_t)Orthorhombic_I];
259 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_ortho_F = lattice_result[(size_t)Orthorhombic_F];
260 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_tetra_P = lattice_result[(size_t)Tetragonal_P];
261 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_tetra_I = lattice_result[(size_t)Tetragonal_I];
262 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_rhom = lattice_result[(size_t)Rhombohedral];
263 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_hex = lattice_result[(size_t)Hexagonal];
264 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_cubic_P = lattice_result[(size_t)Cubic_P];
265 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_cubic_I = lattice_result[(size_t)Cubic_I];
266 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_cubic_F = lattice_result[(size_t)Cubic_F];
267
268 const Int4 num_tri = lattice_result_tri.size();
269
270 /* 2011.10.19 VIC Tamura */
271 Int4 LOOP_COUNTER = 0;
272
273 #ifdef _OPENMP
274 #pragma omp parallel
275 #endif
276 {
277 vector<VCLatticeFigureOfMeritToCheckSymmetry> latFOM_tray;
278
279 #ifdef _OPENMP
280 #pragma omp for
281 #endif
282 for(Int4 n=0; n<num_tri; n++)
283 {
284 /* 2011.10.19 VIC Tamura */
285 SET_PROGRESS(99*(LOOP_COUNTER++)/num_tri, 66, 30); // critical, but works
286 if(IS_CANSELED()) continue;
287
288 VCLatticeFigureOfMeritToCheckSymmetry& latFOM = lattice_result_tri[n];
289 latFOM.setLabel(n+1);
290
291 // Calculate figures of merit as triclinic
292 const ReducedVCLatticeToCheckBravais RLCB(abc_axis, rh_axis, m_DoesPrudentSymSearch, m_cv2, latFOM.putInitialSellingReducedForm() );
293
294 if( JudgeSymmetry[Monoclinic_B] )
295 {
296 putCentringTypes(RLCB, latFOM, BravaisType( put_monoclinic_b_type(abc_axis) ), latFOM_tray);
297 #ifdef _OPENMP
298 #pragma omp critical(monoB)
299 #endif
300 {
301 lattice_result_mono_B.insert(lattice_result_mono_B.end(), latFOM_tray.begin(), latFOM_tray.end());
302 }
303 }
304 if( JudgeSymmetry[Orthorhombic_I] )
305 {
306 putCentringTypes(RLCB, latFOM, BravaisType( pair<eCentringType, ePointGroup>(Inner, D2h) ), latFOM_tray);
307 #ifdef _OPENMP
308 #pragma omp critical(orthoI)
309 #endif
310 {
311 lattice_result_ortho_I.insert(lattice_result_ortho_I.end(), latFOM_tray.begin(), latFOM_tray.end());
312 }
313 }
314 if( JudgeSymmetry[Orthorhombic_F] )
315 {
316 putCentringTypes(RLCB, latFOM, BravaisType( pair<eCentringType, ePointGroup>(Face, D2h) ), latFOM_tray);
317 #ifdef _OPENMP
318 #pragma omp critical(orthoF)
319 #endif
320 {
321 lattice_result_ortho_F.insert(lattice_result_ortho_F.end(), latFOM_tray.begin(), latFOM_tray.end());
322 }
323 }
324 if( JudgeSymmetry[Rhombohedral] )
325 {
326 putCentringTypes(RLCB, latFOM, BravaisType( put_rhombohedral_type(rh_axis) ), latFOM_tray);
327 #ifdef _OPENMP
328 #pragma omp critical(rhom)
329 #endif
330 {
331 lattice_result_rhom.insert(lattice_result_rhom.end(), latFOM_tray.begin(), latFOM_tray.end());
332 }
333 }
334
335 if( JudgeSymmetry[Monoclinic_P] )
336 {
337 latFOM.putLatticesOfHigherSymmetry(put_monoclinic_p_type(abc_axis), m_cv2, latFOM_tray);
338 #ifdef _OPENMP
339 #pragma omp critical(monoP)
340 #endif
341 {
342 lattice_result_mono_P.insert(lattice_result_mono_P.end(), latFOM_tray.begin(), latFOM_tray.end());
343 }
344 }
345 if( JudgeSymmetry[Orthorhombic_P] )
346 {
347 latFOM.putLatticesOfHigherSymmetry(D2h, m_cv2, latFOM_tray);
348 #ifdef _OPENMP
349 #pragma omp critical (ortho_P)
350 #endif
351 {
352 lattice_result_ortho_P.insert(lattice_result_ortho_P.end(), latFOM_tray.begin(), latFOM_tray.end());
353 }
354 }
355 }
356 }
357
358 sort( lattice_result_tri.begin(), lattice_result_tri.end(), VCLatticeFigureOfMeritToCheckSymmetry::cmpFOMdeWolff );
359 if( MAX_SIZE < (Int4)lattice_result_tri.size() )
360 {
361 lattice_result_tri.erase( lattice_result_tri.begin() + MAX_SIZE, lattice_result_tri.end() );
362 }
363
364 const size_t num_tri2 = lattice_result_tri.size();
365 for(size_t n=0; n<num_tri2; n++)
366 {
367 lattice_result_tri[n].setLabel(n+1);
368 }
369
370 ZLOG_INFO( "The program has selected " + num2str<Int4>( lattice_result_tri.size() )
371 + " triclinic solutions using the Wu figure of merit.\n\n" );
372
373
374 /* 2011.10.19 VIC Tamura */
375 CHECK_INTERRUPTION();
376
377 /* 2011.10.19 VIC Tamura */
378 LOOP_COUNTER = 0;
379 Int4 SUM = 0;
380 for(Int4 i=0; i<NUM_LS; i++) { SUM += lattice_result[i].size(); }
381
382 for(Int4 i=1; i<NUM_LS; i++)
383 {
384 if( !JudgeSymmetry[i] ) continue;
385 // sort( lattice_result[i].begin(), lattice_result[i].end() );
386
387 const Int4 num_lattice = lattice_result[i].size();
388
389 #ifdef _OPENMP
390 #pragma omp parallel
391 #endif
392 {
393 vector< VecDat3<Int4> > closest_hkl_tray;
394 vector<bool> is_cal_Q_observed_tray;
395 vector<VCLatticeFigureOfMeritToCheckSymmetry> latFOM_tray;
396
397 #ifdef _OPENMP
398 #pragma omp for
399 #endif
400 for(Int4 index=0; index<num_lattice; index++)
401 {
402 /* 2011.10.19 VIC Tamura */
403 SET_PROGRESS(99+1*(LOOP_COUNTER++)/SUM, 66, 30); // critical, but works
404 if(IS_CANSELED()) continue;
405
406 VCLatticeFigureOfMeritToCheckSymmetry& latFOM0 = lattice_result[i][index];
407 latFOM0.setLabel(index+1);
408
409 latFOM0.setFigureOfMerit(m_num_ref_figure_of_merit,
410 VCData::putPeakQData(),
411 closest_hkl_tray, is_cal_Q_observed_tray);
412
413 LatticeFigureOfMeritZeroShift latFOM2 = latFOM0.putLatticeFigureOfMerit();
414 pair<bool, ZErrorMessage> ans = latFOM2.fitLatticeParameterLinear(VCData::putPeakQData(),
415 closest_hkl_tray, is_cal_Q_observed_tray, false);
416 if( ans.first )
417 {
418 assert( latFOM0.putLatticeFigureOfMerit().putFiguresOfMerit().putNumberOfReflectionsForFigureOfMerit() > 0 );
419 latFOM2.setFigureOfMerit(latFOM0.putLatticeFigureOfMerit().putFiguresOfMerit().putNumberOfReflectionsForFigureOfMerit(),
420 VCData::putPeakQData(),
421 closest_hkl_tray, is_cal_Q_observed_tray);
422 if( LatticeFigureOfMerit::cmpFOMWu( latFOM2, latFOM0.putLatticeFigureOfMerit() ) )
423 {
424 latFOM0.setLatticeFigureOfMerit(latFOM2);
425 }
426 }
427
428 const LatticeFigureOfMerit::SetOfFigureOfMerit& setFOM = latFOM0.putLatticeFigureOfMerit().putFiguresOfMerit();
429 if( setFOM.putFigureOfMeritWu() < MIN_NormM ) continue;
430 if( setFOM.putReversedFigureOfMerit() < MIN_RevM ) continue;
431
432 if( eBravaisType(i) == Monoclinic_P )
433 {
434 if( JudgeSymmetry[Hexagonal] )
435 {
436 latFOM0.putLatticesOfHigherSymmetry(D6h, m_cv2, latFOM_tray);
437 #ifdef _OPENMP
438 #pragma omp critical (hex)
439 #endif
440 {
441 lattice_result_hex.insert(lattice_result_hex.end(), latFOM_tray.begin(), latFOM_tray.end());
442 }
443 }
444 }
445 else if( eBravaisType(i) == Monoclinic_B )
446 {
447 if( JudgeSymmetry[Orthorhombic_C] )
448 {
449 latFOM0.putLatticesOfHigherSymmetry(D2h, m_cv2, latFOM_tray);
450 #ifdef _OPENMP
451 #pragma omp critical (ortho_B)
452 #endif
453 {
454 lattice_result_ortho_B.insert(lattice_result_ortho_B.end(), latFOM_tray.begin(), latFOM_tray.end());
455 }
456 }
457 }
458 else if( eBravaisType(i) == Orthorhombic_P )
459 {
460 if( JudgeSymmetry[Tetragonal_P] )
461 {
462 latFOM0.putLatticesOfHigherSymmetry(D4h_Z, m_cv2, latFOM_tray);
463 #ifdef _OPENMP
464 #pragma omp critical (tetra_P)
465 #endif
466 {
467 lattice_result_tetra_P.insert(lattice_result_tetra_P.end(), latFOM_tray.begin(), latFOM_tray.end());
468 }
469 }
470 if( JudgeSymmetry[Cubic_P] )
471 {
472 latFOM0.putLatticesOfHigherSymmetry(Oh, m_cv2, latFOM_tray);
473 #ifdef _OPENMP
474 #pragma omp critical (cubic_P)
475 #endif
476 {
477 lattice_result_cubic_P.insert(lattice_result_cubic_P.end(), latFOM_tray.begin(), latFOM_tray.end());
478 }
479 }
480 }
481 else if( eBravaisType(i) == Orthorhombic_I )
482 {
483 if( JudgeSymmetry[Tetragonal_I] )
484 {
485 latFOM0.putLatticesOfHigherSymmetry(D4h_Z, m_cv2, latFOM_tray);
486 #ifdef _OPENMP
487 #pragma omp critical (tetra_I)
488 #endif
489 {
490 lattice_result_tetra_I.insert(lattice_result_tetra_I.end(), latFOM_tray.begin(), latFOM_tray.end());
491 }
492 }
493 if( JudgeSymmetry[Cubic_I] )
494 {
495 latFOM0.putLatticesOfHigherSymmetry(Oh, m_cv2, latFOM_tray);
496 #ifdef _OPENMP
497 #pragma omp critical (cubic_I)
498 #endif
499 {
500 lattice_result_cubic_I.insert(lattice_result_cubic_I.end(), latFOM_tray.begin(), latFOM_tray.end());
501 }
502 }
503 }
504 else if( eBravaisType(i) == Orthorhombic_F )
505 {
506 if( JudgeSymmetry[Cubic_F] )
507 {
508 latFOM0.putLatticesOfHigherSymmetry(Oh, m_cv2, latFOM_tray);
509 #ifdef _OPENMP
510 #pragma omp critical (cubic_F)
511 #endif
512 {
513 lattice_result_cubic_F.insert(lattice_result_cubic_F.end(), latFOM_tray.begin(), latFOM_tray.end());
514 }
515 }
516 }
517 }
518 }
519 /* 2011.10.19 VIC Tamura */
520 CHECK_INTERRUPTION();
521
522 sort( lattice_result[i].begin(), lattice_result[i].end(), VCLatticeFigureOfMeritToCheckSymmetry::cmpFOMdeWolff );
523 if( MAX_SIZE < (Int4)lattice_result[i].size() )
524 {
525 lattice_result[i].erase( lattice_result[i].begin() + MAX_SIZE, lattice_result[i].end() );
526 }
527
528 const size_t num_lattice2 = lattice_result[i].size();
529 for(size_t n=0; n<num_lattice2; n++)
530 {
531 lattice_result[i][n].setLabel(n+1);
532 }
533
534 ZLOG_INFO( "(" + num2str( i+1 ) + ") The number of candidates for " + put_bravais_type_name(eBravaisType(i), abc_axis)
535 + " : " + num2str<Int4>( lattice_result[i].size() ) + "\n" );
536 }
537 ZLOG_INFO( "\n" );
538 }
539 catch(bad_alloc& ball)
540 {
541 throw nerror(ball, __FILE__, __LINE__, __FUNCTION__);
542 }
543 }
544
545
546 void SortingLattice::putLatticeCandidatesForEachBravaisTypes(const vector<SymMat43_VCData>& S_super,
547 const Double& MIN_NormM,
548 const Double& MIN_RevM,
549 const Int4& MAX_SIZE,
550 const eABCaxis& abc_axis,
551 const eRHaxis& rh_axis,
552 vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[NUM_LS]) const
553 {
554 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result_tri = lattice_result[(size_t)Triclinic];
555 putLatticeCandidatesForTriclinic(S_super, MIN_NormM, MIN_RevM, lattice_result_tri);
556
557 ZLOG_INFO( "Determination of the Bravais type is being carried out...\n(Solutions with " + putLabel(SCWuM) + " less than " + num2str(MIN_NormM)
558 + " or " + putLabel(SCRevM) + " less than " + num2str(MIN_RevM)
559 + " are automatically removed).\n" );
560 ZLOG_INFO( "All the unit-cell parameters are being optimized by linear least squares...\n" );
561
562 //ZLOG_INFO( "The program has removed " + num2str<Int4>( S_super.size() - lattice_result_tri.size() )
563 // + " triclinic solutions because their " + putLabel(SCWuM) + " is less than " + num2str(MIN_NormM)
564 // + " or their " + putLabel(SCRevM) + " is less than " + num2str(MIN_RevM) + ".\n\n" );
565 //ZLOG_INFO( "Determination of the Bravais type is being carried out...\n" );
566 putLatticeCandidatesForEachBravaisTypes(MIN_NormM, MIN_RevM, MAX_SIZE, abc_axis, rh_axis, lattice_result);
567 }
568
569
570 void SortingLattice::setNumberOfNeighbors(const eABCaxis& baxis_flag,
571 bool (*CmpFunc)(const VCLatticeFigureOfMeritToCheckSymmetry&, const VCLatticeFigureOfMeritToCheckSymmetry&),
572 vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[NUM_LS]) const
573 {
574
575 #ifdef _OPENMP
576 #pragma omp for
577 #endif
578 for(Int4 i=0; i<NUM_LS; i++)
579 {
580 if( !OutputSymmetry[(size_t)i] ) continue;
581
582 stable_sort( lattice_result[i].begin(), lattice_result[i].end() ); // Sort by the unit-cell volume.
583 for(vector<VCLatticeFigureOfMeritToCheckSymmetry>::iterator it=lattice_result[i].begin(); it<lattice_result[i].end(); it++)
584 {
585 it->putNumberOfLatticesInNeighborhood() = 0;
586 }
587 }
588
589 const Double coef_lower = 1.0 - m_resol*3.0;
590 const Double coef_upper = 1.0 + m_resol*3.0;
591
592 Vec_INT index_tray(put_number_of_bravais_types(), 0);
593
594 /* 2011.10.19 VIC Tamura */
595 Int4 SUM=0, LOOP_COUNTER=0;
596 for(Int4 i=0; i<NUM_LS; i++ ) { SUM += lattice_result[(size_t)i].size(); }
597
598 #ifdef _OPENMP
599 #pragma omp for
600 #endif
601 for(Int4 i=0; i<NUM_LS; i++)
602 {
603 if( !OutputSymmetry[(size_t)i] ) continue;
604
605 const size_t num_lattice = lattice_result[i].size();
606
607 for(size_t index=0; index<num_lattice; index++)
608 {
609 /* 2011.10.19 VIC Tamura */
610 SET_PROGRESS(100*(LOOP_COUNTER++)/SUM, 97, 1); // critical, but works
611 if(IS_CANSELED()) continue;
612
613 VCLatticeFigureOfMeritToCheckSymmetry& latFOM0 = lattice_result[i][index];
614 const LatticeFigureOfMerit& latFOM0_prim = latFOM0.putLatticeFigureOfMerit();
615 if( latFOM0.putNumberOfLatticesInNeighborhood() < 0 ) continue;
616
617 const Double& detS = latFOM0_prim.putDeterminantOfGramMatrix();
618 const size_t ibegin = distance( lattice_result[i].begin(), lower_bound( lattice_result[i].begin(), lattice_result[i].end(), detS*coef_lower ) );
619 const size_t iend = distance( lattice_result[i].begin(), upper_bound( lattice_result[i].begin()+ibegin, lattice_result[i].end(), detS*coef_upper ) );
620
621 Int4 count=0;
622 if( i == (size_t)Triclinic )
623 {
624 const ReducedLatticeToCheckEquiv RLCS(m_resol, latFOM0_prim.putSellingReducedForm());
625 for(size_t index2=ibegin; index2<iend; index2++)
626 {
627 if( index2 == index ) continue;
628
629 VCLatticeFigureOfMeritToCheckSymmetry& latFOM2 = lattice_result[i][index2];
630 const LatticeFigureOfMerit& latFOM2_prim = latFOM2.putLatticeFigureOfMerit();
631
632 // lattice_result_tri[index2] equals trans_mat * RLCB.m_S_super_obtuse * Transpose(trans_mat)
633 if( RLCS.equiv( latFOM2_prim.putSellingReducedForm() ) )
634 {
635 // Compare the figures of merit.
636 if( CmpFunc( latFOM2, latFOM0 ) )
637 {
638 if( latFOM2.putNumberOfLatticesInNeighborhood() >= 0 )
639 {
640 latFOM2.putNumberOfLatticesInNeighborhood() += 1;
641 count = -1;
642 break;
643 }
644 }
645 else
646 {
647 count++;
648 latFOM2.putNumberOfLatticesInNeighborhood() = -1;
649 }
650 }
651 }
652 }
653 else
654 {
655 for(size_t index2=ibegin; index2<iend; index2++)
656 {
657 if( index2 == index ) continue;
658
659 VCLatticeFigureOfMeritToCheckSymmetry& latFOM2 = lattice_result[i][index2];
660 const LatticeFigureOfMerit& latFOM2_prim = latFOM2.putLatticeFigureOfMerit();
661
662 // *it2 equals trans_mat * RLCS.m_S_super_obtuse * Transpose(trans_mat)
663 if( check_equiv_m( latFOM0_prim.putInverseOfBuergerReducedForm(), latFOM2_prim.putInverseOfBuergerReducedForm(), m_resol ) )
664 {
665 // Compare the figures of merit.
666 if( CmpFunc( latFOM2, latFOM0 ) )
667 {
668 if( latFOM2.putNumberOfLatticesInNeighborhood() >= 0 )
669 {
670 latFOM2.putNumberOfLatticesInNeighborhood() += 1;
671 count = -1;
672 break;
673 }
674 }
675 else
676 {
677 count++;
678 latFOM2.putNumberOfLatticesInNeighborhood() = -1;
679 }
680 }
681 }
682 }
683
684 latFOM0.putNumberOfLatticesInNeighborhood() = count;
685 }
686
687 Int4& index = index_tray[i];
688 index = 0;
689 for(vector<VCLatticeFigureOfMeritToCheckSymmetry>::const_iterator it=lattice_result[i].begin(); it<lattice_result[i].end(); it++)
690 {
691 if( it->putNumberOfLatticesInNeighborhood() >= 0 ) index++;
692 }
693 }
694
695 /* 2011.10.19 VIC Tamura */
696 CHECK_INTERRUPTION();
697
698 for(Int4 i=0; i<NUM_LS; i++)
699 {
700 if( !OutputSymmetry[(size_t)i] ) continue;
701 ZLOG_INFO( "(" + num2str( i+1 ) + ") The number of candidates for " + put_bravais_type_name(eBravaisType(i), baxis_flag)
702 + " : " + num2str( lattice_result[i].size() ) + " -> " + num2str( index_tray[i] ) + "\n" );
703 }
704 ZLOG_INFO( "\n" );
705 }

Back to OSDN">Back to OSDN
ViewVC Help
Powered by ViewVC 1.1.26