Develop and Download Open Source Software

Powder auto-indexing software Conograph

Browse Subversion Repository

/[conograph]/Conograph/trunk/src/p_out_indexing.cc

Contents of /Conograph/trunk/src/p_out_indexing.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: 26335 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 #include <iostream>
28 #include "p_out_indexing.hh"
29 #include "ControlParam.hh"
30 #include "PeakPosData.hh"
31 #include "lattice_symmetry/VCLatticeFigureOfMeritToCheckSymmetry.hh"
32 #include "lattice_symmetry/LatticeFigureOfMeritToDisplay.hh"
33 #include "utility_func/covar_matrix.hh"
34 #include "utility_func/zmath.hh"
35 #include "utility_lattice_reduction/matrix_NbyN.hh"
36 #include "utility_data_structure/Bud.hh"
37 #include "utility_data_structure/Node3.hh"
38 #include "utility_data_structure/qdata.hh"
39 #include "utility_data_structure/SymMatWCovar.hh"
40 #include "utility_data_structure/VecDat3.hh"
41 #include "lattice_symmetry/OutputInfo.hh"
42
43
44 void printQdata(const vector<QData>& qdata,
45 ostream* os)
46 {
47 os->setf(ios::right);
48 os->setf(ios::uppercase);
49 os->setf(ios::showpoint);
50 os->setf(ios::scientific);
51
52 *os << "IGOR" << endl;
53 *os << "WAVES/O index, Q, sqrt_root_of_variance_of_Q" << endl;
54 *os << "BEGIN" << endl;
55
56 const Int4 isize = qdata.size();
57
58 for (Int4 k=0; k<isize; k++)
59 {
60 os->width(5);
61 *os << k+1;
62
63 os->width(15);
64 *os << qdata[k].q;
65
66 os->width(15);
67 *os << sqrt(qdata[k].q_var);
68
69 *os << endl;
70 }
71 *os << "END" << endl;
72 }
73
74
75
76 void print_bud_data(const vector<Bud>& bud_basket,
77 const Vec_BOOL& root_flag,
78 const vector< vector<Int4> >& leftbr_tray,
79 const vector< vector<Int4> >& rightbr_tray,
80 const vector< vector<Int4> >& root_leftbr,
81 const vector< vector<Int4> >& root_rightbr,
82 const string& fname)
83 {
84 ofstream ofs(fname.c_str());
85 ostream *os;
86 os = &ofs;
87
88 *os << "IGOR" << endl;
89 *os << "WAVES/O ";
90
91 *os << "Number, K1, K2, K3, K4, (Q1+Q2)*2-Q3-Q4, Err[(Q1+Q2)*2-Q3-Q4]";
92 *os << "root?, super_base?, left_branch, right_branch, left_root, right_root, tree_index" << endl;
93
94 *os << "BEGIN" << endl;
95 os->setf(ios::right);
96 os->setf(ios::uppercase);
97 os->setf(ios::showpoint);
98
99 Int4 k=0;
100 for(vector<Bud>::const_iterator it=bud_basket.begin(); it!=bud_basket.end(); it++, k++)
101 {
102 os->unsetf(ios::scientific);
103 os->width(5);
104 *os << k+1;
105
106 os->width(5);
107 if( it->iK1() < 0 ) *os << -1;
108 else *os << it->iK1() + 1;
109
110 os->width(5);
111 if( it->iK2() < 0 ) *os << -1;
112 else *os << it->iK2() + 1;
113
114 os->width(5);
115 if( it->iK3() < 0 ) *os << -1;
116 else *os << it->iK3() + 1;
117
118 os->width(5);
119 if( it->iK4() < 0 ) *os << -1;
120 else *os << it->iK4() + 1;
121
122 const VCData diff = (it->Q1()+it->Q2())*2-(it->Q3()+it->Q4());
123
124 os->setf(ios::scientific);
125 os->width(15);
126 *os << diff.Value();
127
128 os->width(15);
129 *os << sqrt( diff.Variance() );
130
131 os->precision(6);
132 os->width(5);
133 if( root_flag[k] ) *os << 1;
134 else *os << 0;
135
136 os->width(5);
137 if( it->IsSuperBasisObtuse() ) *os << 1;
138 else *os << 0;
139
140 os->width(10);
141 if( leftbr_tray[k].empty() ) *os << "Emp";
142 else
143 {
144 *os << leftbr_tray[k][0] + 1;
145 for(UInt4 l=1; l<leftbr_tray[k].size(); l++) *os << ":" << leftbr_tray[k][l] + 1;
146 }
147
148 os->width(10);
149 if( rightbr_tray[k].empty() ) *os << "Emp";
150 else
151 {
152 *os << rightbr_tray[k][0] + 1;
153 for(UInt4 l=1; l<rightbr_tray[k].size(); l++) *os << ":" << rightbr_tray[k][l] + 1;
154 }
155
156 os->width(10);
157 if( root_leftbr[k].empty() ) *os << "Emp";
158 else
159 {
160 *os << root_leftbr[k][0] + 1;
161 for(UInt4 l=1; l<root_leftbr[k].size(); l++) *os << ":" << root_leftbr[k][l] + 1;
162 }
163
164 os->width(10);
165 if( root_rightbr[k].empty() ) *os << "Emp";
166 else
167 {
168 *os << root_rightbr[k][0] + 1;
169 for(UInt4 l=1; l<root_rightbr[k].size(); l++) *os << ":" << root_rightbr[k][l] + 1;
170 }
171
172 *os << endl;
173 }
174 *os << "END" << endl;
175
176 ofs.close();
177 }
178
179
180
181
182 void print_bud_data(const vector<Bud>& bud_basket, const string& fname)
183 {
184 ofstream ofs(fname.c_str());
185 ostream *os;
186 os = &ofs;
187
188 *os << "IGOR" << endl;
189 *os << "WAVES/O ";
190
191 *os << "Number, K1, K2, K3, K4, (Q1+Q2)*2-Q3-Q4, Err[(Q1+Q2)*2-Q3-Q4]";
192 *os << "root?, super_base?, left_branch, right_branch, left_root, right_root, tree_index" << endl;
193
194 *os << "BEGIN" << endl;
195 os->setf(ios::right);
196 os->setf(ios::uppercase);
197 os->setf(ios::showpoint);
198
199 Int4 k=0;
200 for(vector<Bud>::const_iterator it=bud_basket.begin(); it!=bud_basket.end(); it++, k++)
201 {
202 os->unsetf(ios::scientific);
203 os->width(5);
204 *os << k+1;
205
206 os->width(5);
207 if( it->iK1() < 0 ) *os << -1;
208 else *os << it->iK1() + 1;
209
210 os->width(5);
211 if( it->iK2() < 0 ) *os << -1;
212 else *os << it->iK2() + 1;
213
214 os->width(5);
215 if( it->iK3() < 0 ) *os << -1;
216 else *os << it->iK3() + 1;
217
218 os->width(5);
219 if( it->iK4() < 0 ) *os << -1;
220 else *os << it->iK4() + 1;
221
222 const VCData diff = (it->Q1()+it->Q2())*2-(it->Q3()+it->Q4());
223
224 os->setf(ios::scientific);
225 os->width(15);
226 *os << diff.Value();
227
228 os->width(15);
229 *os << sqrt( diff.Variance() );
230
231 os->width(5);
232 if( it->IsSuperBasisObtuse() ) *os << 1;
233 else *os << 0;
234
235 *os << endl;
236 }
237 *os << "END" << endl;
238
239 ofs.close();
240 }
241
242
243
244
245 void print_node_data(const vector<Node3>& node_basket,
246 const string& fname)
247 {
248 ofstream ofs(fname.c_str());
249 ostream *os;
250 os = &ofs;
251
252 ofs << "** Node, K1, K2, K3, K1-K2, K1-K3, K2-K3, Q1, Q2, Q3, Q12, Q13, Q23, UnitcellVolume" << endl;
253
254 os->setf(ios::scientific);
255 os->precision(4);
256
257 Int4 iK1, iK2, iK3;
258 Int4 iK12, iK13, iK23;
259 Double detS;
260 for(UInt4 l=0; l<node_basket.size(); l++)
261 {
262 os->width(5);
263 *os << l + 1;
264
265 const Node3& nodex = node_basket[l];
266 iK1 = nodex.iK(0);
267 iK2 = nodex.iK(1);
268 iK3 = nodex.iK(2);
269
270 iK12 = nodex.iSumK(0,1);
271 iK13 = nodex.iSumK(0,2);
272 iK23 = nodex.iSumK(1,2);
273
274 os->width(8);
275 if( iK1 >= 0 ) *os << iK1 + 1;
276 else *os << -1;
277
278 os->width(5);
279 if( iK2 >= 0 ) *os << iK2 + 1;
280 else *os << -1;
281
282 os->width(5);
283 if( iK3 >= 0 ) *os << iK3 + 1;
284 else *os << -1;
285
286 os->width(8);
287 if( iK12 >= 0 ) *os << iK12 + 1;
288 else *os << -1;
289
290 os->width(5);
291 if( iK13 >= 0 ) *os << iK13 + 1;
292 else *os << -1;
293
294 os->width(5);
295 if( iK23 >= 0 ) *os << iK23 + 1;
296 else *os << -1;
297
298
299 os->width(15);
300 *os << nodex.K(0).Value();
301
302 os->width(15);
303 *os << nodex.K(1).Value();
304
305 os->width(15);
306 *os << nodex.K(2).Value();
307
308 os->width(15);
309 *os << nodex.SumK(0,1).Value();
310
311 os->width(15);
312 *os << nodex.SumK(0,2).Value();
313
314 os->width(15);
315 *os << nodex.SumK(1,2).Value();
316
317 detS = nodex.putDeterminantS();
318
319 os->width(20);
320 *os << 1.0/sqrt(detS);
321
322 *os << endl;
323 }
324 *os << endl;
325 }
326
327
328
329
330 void printHKLdata(const vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[],
331 const OutputInfo outinfo[],
332 const eSortCriterion& sort_criterion,
333 const ControlParam& cdata,
334 const PeakPosData& pdata,
335 const string& fname)
336 {
337 static const Int4 NUM_LS = put_number_of_bravais_types();
338
339 ofstream ofs(fname.c_str());
340 ostream *os;
341 os = &ofs;
342
343 os->setf(ios::scientific);
344 os->precision(4);
345
346 pair< vector<Int4>::const_iterator, vector<Int4>::const_iterator> it_pair;
347 SymMatWCovar dbl_S(3);
348 VecDat3<Double> length_axis, angle_axis;
349 SymMat<VCData> S_red2(3);
350 set< SymMat<VCData> > S_tray;
351
352 Int4 label_start=0;
353 os->width(label_start);
354 *os << "" << "<?xml version=\"1.0\" encoding=\"UTF-8\" ?>\n";
355
356 os->width(label_start);
357 *os << "" << "<ZCodeParameters>\n";
358 label_start++;
359
360 os->width(label_start);
361 *os << "" << "<ConographOutput>\n";
362 label_start++;
363
364 os->width(label_start);
365 *os << "" << "<!-- Information on the best " + putLabel(sort_criterion) + " solution for each Bravais type.\n";
366 os->width(label_start+6);
367 *os << "" << "TNB : total number of solutions of the Bravais types,\n";
368 os->width(label_start+6);
369 *os << "" << putLabel(SCM) << " : de Wolff figure of merit,\n";
370 os->width(label_start+6);
371 *os << "" << putLabel(SCWuM) << " : Wu figure of merit,\n";
372 // os->width(label_start+6);
373 // *os << "" << putLabel(SCNormM) << " : normalized de Wolff figure of merit,\n";
374 os->width(label_start+6);
375 *os << "" << putLabel(SCRevM) << " : reversed de Wolff figure of merit,\n";
376 os->width(label_start+6);
377 *os << "" << putLabel(SCSymM) << " : symmetric de Wolff figure of merit,\n";
378 os->width(label_start+6);
379 *os << "" << putLabel(SCNN) << " : number of lattices in the neighborhood,\n";
380 os->width(label_start+6);
381 *os << "" << "VOL : unit-cell volume.\n";
382 os->width(label_start);
383 *os << "" << "Bravais Lattice : TNB";
384 for(Int4 i=0; i<putNumberOfSortCriterion(); i++)
385 {
386 *os << ", " << putLabel(eSortCriterion(i));
387 }
388 *os << ", VOL\n";
389
390 for(Int4 k=NUM_LS-1; k>=0; k--)
391 {
392 os->width(17);
393 *os << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << " : ";
394 os->width(7);
395 *os << outinfo[k].putCandidateNumToOutput();
396
397 const Int4 index = outinfo[k].putIndex( outinfo[k].putLatticeLabelSelectedFromListToOutput(sort_criterion) );
398 if( index >= 0 )
399 {
400 const LatticeFigureOfMerit& lat_fom = lattice_result[k][index].putLatticeFigureOfMerit();
401 const LatticeFigureOfMerit::SetOfFigureOfMerit& set_fom = lat_fom.putFiguresOfMerit();
402
403 os->width(15);
404 *os << set_fom.putFigureOfMeritWolff();
405
406 os->width(15);
407 *os << set_fom.putFigureOfMeritWu();
408
409 // os->width(15);
410 // *os << set_fom.putNormalizedFigureOfMeritWolff();
411
412 os->width(15);
413 *os << set_fom.putReversedFigureOfMerit();
414
415 os->width(15);
416 *os << set_fom.putSymmetricFigureOfMerit();
417
418 os->width(5);
419 *os << lattice_result[k][index].putNumberOfLatticesInNeighborhood();
420
421 os->width(15);
422 *os << lat_fom.putLatticeVolume();
423
424 // os->width(15);
425 // *os << set_fom.putFigureOfMeritWolff_Original();
426 //
427 // os->width(15);
428 // *os << set_fom.putFigureOfMeritWu_Original();
429 }
430 *os << endl;
431 }
432 os->width(label_start);
433 *os << "" << "-->\n\n";
434 os->width(label_start);
435 *os << "" << "<!-- Labels of the solution with the best figure of merit.\n";
436 for(Int4 k=NUM_LS-1; k>=0; k--)
437 {
438 os->width(18);
439 *os << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << ", Best Score : Lattice constants, label.\n";
440 // if( lattice_result[k].empty() )
441 // {
442 // *os << "\n";
443 // continue;
444 // }
445 for(Int4 i=0; i<putNumberOfSortCriterion(); i++)
446 {
447 os->width(15);
448 *os << putLabel(eSortCriterion(i)) << " = ";
449
450 const Int4 index = outinfo[k].putIndex( outinfo[k].putLatticeLabelSelectedAmongAll(eSortCriterion(i)) );
451 os->width(12);
452 if( index >= 0 )
453 {
454 const LatticeFigureOfMerit& lat_fom = lattice_result[k][index].putLatticeFigureOfMerit();
455 const LatticeFigureOfMerit::SetOfFigureOfMerit& set_fom = lat_fom.putFiguresOfMerit();
456
457 if( eSortCriterion(i) == SCM )
458 {
459 *os << set_fom.putFigureOfMeritWolff();
460 }
461 else if( eSortCriterion(i) == SCWuM )
462 {
463 *os << set_fom.putFigureOfMeritWu();
464 }
465 // if( eSortCriterion(i) == SCNormM )
466 // {
467 // *os << set_fom.putNormalizedFigureOfMeritWolff();
468 // }
469 else if( eSortCriterion(i) == SCRevM )
470 {
471 *os << set_fom.putReversedFigureOfMerit();
472 }
473 else if( eSortCriterion(i) == SCSymM )
474 {
475 *os << set_fom.putSymmetricFigureOfMerit();
476 }
477 else if( eSortCriterion(i) == SCNN )
478 {
479 *os << lattice_result[k][index].putNumberOfLatticesInNeighborhood();
480 }
481 // else if( eSortCriterion(i) == SCM_org )
482 // {
483 // *os << set_fom.putFigureOfMeritWolff_Original();
484 // }
485 // else if( eSortCriterion(i) == SCWuM_org )
486 // {
487 // *os << set_fom.putFigureOfMeritWu_Original();
488 // }
489
490 *os << " : ";
491 lat_fom.putReducedLatticeConstantsDegree(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), length_axis, angle_axis);
492 os->width(14);
493 *os << length_axis[0];
494 os->width(14);
495 *os << length_axis[1];
496 os->width(14);
497 *os << length_axis[2];
498 os->width(14);
499 *os << angle_axis[0];
500 os->width(14);
501 *os << angle_axis[1];
502 os->width(14);
503 *os << angle_axis[2];
504 os->width(10);
505 *os << lattice_result[k][index].putStringLabel();
506 *os << endl;
507 }
508 else
509 {
510 *os << "- -" << " : - -\n";
511 }
512 }
513 *os << endl;
514 }
515 os->width(label_start);
516 *os << "" << "-->\n\n";
517
518 pair<eBravaisType, lattice_label> selected_lattice_info = put_selected_lattice_from_all(lattice_result, outinfo, sort_criterion);
519 const Int4 selected_lattice_index = outinfo[selected_lattice_info.first].putIndex(selected_lattice_info.second);
520
521 if( selected_lattice_index >= 0 )
522 {
523 const VCLatticeFigureOfMeritToCheckSymmetry& selected_lattice = lattice_result[selected_lattice_info.first][selected_lattice_index];
524
525 os->width(label_start);
526 *os << "" << "<!-- Information on the selected candidates.-->\n";
527
528 os->width(label_start);
529 *os << "" << "<SelectedLatticeCandidate number=\"" << selected_lattice.putStringLabel() << "\">\n";
530 label_start++;
531
532 selected_lattice.printLatticeInformation(lattice_result, outinfo,
533 cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
534 cdata.putResolution(), label_start, os);
535
536 label_start--;
537 os->width(label_start);
538 *os << "" << "</SelectedLatticeCandidate>\n\n";
539 }
540
541 for(Int4 k=NUM_LS-1; k>=0; k--)
542 {
543 *os << "\n";
544 os->width(label_start);
545 *os << "" << "<!-- Candidates for " << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << " -->\n\n";
546
547 const Int4 num_topo = lattice_result[k].size();
548
549 for(Int4 n=0; n<num_topo; n++)
550 {
551 const VCLatticeFigureOfMeritToCheckSymmetry& ans = lattice_result[k][n];
552 if( !outinfo[k].IsOutput(ans.putLatticeLabel()) ) continue;
553
554 os->width(label_start);
555 *os << "" << "<LatticeCandidate number=\"" << ans.putStringLabel() << "\">\n";
556 label_start++;
557 ans.printLatticeInformation(lattice_result, outinfo, cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
558 0.06, label_start, os);
559
560 label_start--;
561 os->width(label_start);
562 *os << "" << "</LatticeCandidate>\n\n";
563 }
564 *os << endl;
565 }
566
567 label_start--;
568 os->width(label_start);
569 *os << "" << "</ConographOutput>\n";
570
571 label_start--;
572 os->width(label_start);
573 *os << "" << "</ZCodeParameters>\n";
574 }
575
576
577
578 void printHKLdata(const vector<LatticeFigureOfMeritToDisplay> lattice_result,
579 const ControlParam& cdata,
580 const PeakPosData& pdata,
581 const string& fname)
582 {
583 ofstream ofs(fname.c_str());
584 ostream *os;
585 os = &ofs;
586
587 os->setf(ios::scientific);
588 os->precision(4);
589
590 pair< vector<Int4>::const_iterator, vector<Int4>::const_iterator> it_pair;
591 SymMatWCovar dbl_S(3);
592 VecDat3<Double> length_axis, angle_axis;
593 SymMat<VCData> S_red2(3);
594 set< SymMat<VCData> > S_tray;
595
596 Int4 label_start=0;
597 os->width(label_start);
598 *os << "" << "<ZCodeParameters>\n";
599 label_start++;
600
601 os->width(label_start);
602 *os << "" << "<ConographOutput>\n";
603 label_start++;
604
605 const Int4 num_topo = lattice_result.size();
606
607 for(Int4 n=0; n<num_topo; n++)
608 {
609 const LatticeFigureOfMeritZeroShift& ans = lattice_result[n].putLatticeFigureOfMerit();
610
611 os->width(label_start);
612 *os << "" << "<LatticeCandidate>\n";
613 label_start++;
614 ans.printLatticeInformation(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
615 cdata.putResolution(), label_start, os);
616 lattice_result[n].printIndexingResult(cdata, pdata, label_start, os);
617
618 label_start--;
619 os->width(label_start);
620 *os << "" << "</LatticeCandidate>\n\n";
621 }
622 *os << endl;
623
624 label_start--;
625 os->width(label_start);
626 *os << "" << "</ConographOutput>\n";
627
628 label_start--;
629 os->width(label_start);
630 *os << "" << "</ZCodeParameters>\n";
631 }
632
633
634
635
636 void printPeakPosition(
637 const ControlParam& cdata,
638 const PeakPosData& pdata,
639 const LatticeFigureOfMeritToDisplay& latfit,
640 const vector<LatticeFigureOfMerit>& lattices_same_q,
641 const string& fname)
642 {
643 ofstream ofs(fname.c_str());
644 ostream * const os = &ofs;
645
646 *os << "IGOR\n";
647 *os << "WAVES/O ";
648
649 pdata.printData(os);
650
651 *os << "WAVES/O dphase_0, xphase_0, yphase_0, h_0, k_0, l_0\n";
652
653 *os << "BEGIN\n";
654
655 const vector<HKL_Q>& cal_hkl_tray = latfit.putCalMillerIndices();
656 Vec_DP cal_pos_tray;
657 latfit.putCalculatedPeakPosInRange(cdata, cal_pos_tray);
658
659 const Int4 peak_num = cal_hkl_tray.size();
660
661 for (Int4 i=0; i<peak_num; i++)
662 {
663 os->width(15);
664 *os << 1.0 / sqrt( cal_hkl_tray[i].Q() );
665
666 os->width(15);
667 if( cal_pos_tray[i] < 0.0 )
668 {
669 *os << "NAN";
670 }
671 else
672 {
673 *os << cal_pos_tray[i];
674 }
675 os->width(15);
676 *os << 0.0;
677
678 const VecDat3<Int4>& hkl = cal_hkl_tray[i].HKL();
679
680 os->width(5);
681 *os << hkl[0];
682 os->width(5);
683 *os << hkl[1];
684 os->width(5);
685 *os << hkl[2];
686
687 *os << endl;
688 }
689 *os << "END\n\n";
690
691 LatticeFigureOfMeritToDisplay latfit2;
692 latfit2.setPeakShiftParamDegree(latfit.putPeakShiftFunctionType(), latfit.putWaveLength(),
693 latfit.putPeakShiftParamDegree(), pdata);
694
695 const Int4 isize = lattices_same_q.size();
696 for (Int4 j=0, j2=1; j<isize; j++, j2++)
697 {
698 VecDat3<Double> length_axis, angle_axis;
699 lattices_same_q[j].putReducedLatticeConstantsDegree(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), length_axis, angle_axis);
700
701 latfit2.setLatticeConstantsDegree(lattices_same_q[j].enumBravaisType(), cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), length_axis, angle_axis);
702 latfit2.resetMillerIndicesInRange(cdata.putNumberOfReflectionsForFigureOfMerit());
703
704
705 *os << "WAVES/O dphase_" << j2 << ", xphase_" << j2 << ", yphase_" << j2 << ", ";
706 *os << "h_" << j2 << ", ";
707 *os << "k_" << j2 << ", ";
708 *os << "l_" << j2;
709 *os << endl;
710
711 *os << "BEGIN\n";
712
713 const vector<HKL_Q>& cal_hkl_tray = latfit2.putCalMillerIndices();
714 Vec_DP cal_pos_tray;
715 latfit2.putCalculatedPeakPosInRange(cdata, cal_pos_tray);
716
717 const Int4 peak_num = cal_hkl_tray.size();
718
719 for (Int4 i=0; i<peak_num; i++)
720 {
721 os->width(15);
722 *os << 1.0 / sqrt( cal_hkl_tray[i].Q() );
723
724 os->width(15);
725 if( cal_pos_tray[i] < 0.0 )
726 {
727 *os << "NAN";
728 }
729 else
730 {
731 *os << cal_pos_tray[i];
732 }
733 os->width(15);
734 *os << 0.0;
735
736 const VecDat3<Int4>& hkl = cal_hkl_tray[i].HKL();
737
738 os->width(5);
739 *os << hkl[0];
740 os->width(5);
741 *os << hkl[1];
742 os->width(5);
743 *os << hkl[2];
744
745 *os << endl;
746 }
747 *os << "END\n\n";
748 }
749
750 VecDat3<Double> length_axis, angle_axis;
751 const string str_num_ref = num2str( cdata.putNumberOfReflectionsForFigureOfMerit() );
752
753
754 *os << "X Display " << pdata.putYIntColumnTitle() << " vs " << pdata.putXColumnTitle() << endl;
755 *os << "X ModifyGraph mirror(left)=2\n";
756 *os << "X ModifyGraph mirror(bottom)=2\n";
757 *os << "X ModifyGraph rgb(" << pdata.putYIntColumnTitle() << ")=(0,15872,65280)\n";
758
759 *os << "X AppendToGraph height vs peakpos\n";
760 *os << "X ModifyGraph mode(height)=3,marker(height)=17\n";
761 *os << "X ModifyGraph rgb(height)=(0,65280,65280)\n";
762
763 Double offset = 0.0;
764 const Double offset_gap = pdata.putMaxPeakHeightOfFirst20() * (-0.05);
765 for (Int4 j=0; j<isize+1; j++, offset+=offset_gap)
766 {
767 *os << "X AppendToGraph yphase_" << j << " vs xphase_" << j << endl;
768 *os << "X ModifyGraph offset(yphase_" << j << ")={0," << offset << "},mode(yphase_" << j << ")=3,marker(yphase_" << j;
769 *os << ")=10,msize(yphase_" << j << ")=3,mrkThick(yphase_" << j << ")=0.6,rgb(yphase_" << j << ")=(3,52428,1)" << endl;
770 if( j <= 0 ) offset+=offset_gap;
771 }
772
773 *os << "X Legend/C/N=text0/J/A=MC \"";
774 *os << "\\s(yphase_" << 0 << ") "
775 + put_bravais_type_name(latfit.putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
776 + latfit.putLatticeFigureOfMerit().printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
777 for (Int4 j=0, j2=1; j<isize; j++, j2++)
778 {
779 *os << "\\r\\s(yphase_" << j2 << ") "
780 + put_bravais_type_name(lattices_same_q[j].enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
781 + lattices_same_q[j].printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
782 }
783 *os << "\"\nX Legend/C/N=text0/J/A=RT/X=0.00/Y=0.00\n";
784
785 ofs.close();
786 }
787
788
789
790 void printPeakPosition(
791 const ControlParam& cdata,
792 const PeakPosData& pdata,
793 const vector<LatticeFigureOfMeritToDisplay>& latfit_tray,
794 const string& fname)
795 {
796 ofstream ofs(fname.c_str());
797 ostream * const os = &ofs;
798
799 *os << "IGOR\n";
800 *os << "WAVES/O ";
801
802 pdata.printData(os);
803
804 if( latfit_tray.empty() ) return;
805
806 const Int4 isize = latfit_tray.size();
807 for (Int4 j=0, j1=1; j<isize; j++, j1++)
808 {
809 const LatticeFigureOfMeritToDisplay& latfit = latfit_tray[j];
810
811 *os << "WAVES/O dphase_" << j1 << ", xphase_" << j1 << ", yphase_" << j1 << ", ";
812 *os << "h_" << j1 << ", ";
813 *os << "k_" << j1 << ", ";
814 *os << "l_" << j1;
815 *os << endl;
816
817 *os << "BEGIN\n";
818
819 const vector<HKL_Q>& cal_hkl_tray = latfit.putCalMillerIndices();
820 Vec_DP cal_pos_tray;
821 latfit.putCalculatedPeakPosInRange(cdata, cal_pos_tray);
822
823 const Int4 peak_num = cal_hkl_tray.size();
824
825 for (Int4 i=0; i<peak_num; i++)
826 {
827 os->width(15);
828 *os << 1.0 / sqrt( cal_hkl_tray[i].Q() );
829
830 os->width(15);
831 if( cal_pos_tray[i] < 0.0 )
832 {
833 *os << "NAN";
834 }
835 else
836 {
837 *os << cal_pos_tray[i];
838 }
839 os->width(15);
840 *os << 0.0;
841
842 const VecDat3<Int4>& hkl = cal_hkl_tray[i].HKL();
843
844 os->width(5);
845 *os << hkl[0];
846 os->width(5);
847 *os << hkl[1];
848 os->width(5);
849 *os << hkl[2];
850
851 *os << endl;
852 }
853 *os << "END\n\n";
854 }
855
856 VecDat3<Double> length_axis, angle_axis;
857 const string str_num_ref = num2str( cdata.putNumberOfReflectionsForFigureOfMerit() );
858
859
860 *os << "X Display " << pdata.putYIntColumnTitle() << " vs " << pdata.putXColumnTitle() << endl;
861 *os << "X ModifyGraph mirror(left)=2\n";
862 *os << "X ModifyGraph mirror(bottom)=2\n";
863 *os << "X ModifyGraph rgb(" << pdata.putYIntColumnTitle() << ")=(0,15872,65280)\n";
864
865 *os << "X AppendToGraph height vs peakpos\n";
866 *os << "X ModifyGraph mode(height)=3,marker(height)=17\n";
867 *os << "X ModifyGraph rgb(height)=(0,65280,65280)\n";
868
869 Double offset = 0.0;
870 const Double offset_gap = pdata.putMaxPeakHeightOfFirst20() * (-0.05);
871 for (Int4 j=0, j1=1; j<isize; j++, j1++, offset+=offset_gap)
872 {
873 *os << "X AppendToGraph yphase_" << j1 << " vs xphase_" << j1 << endl;
874 *os << "X ModifyGraph offset(yphase_" << j1 << ")={0," << offset << "},mode(yphase_" << j1 << ")=3,marker(yphase_" << j1;
875 *os << ")=10,msize(yphase_" << j1 << ")=3,mrkThick(yphase_" << j1 << ")=0.6,rgb(yphase_" << j1 << ")=(3,52428,1)" << endl;
876 }
877
878 *os << "X Legend/C/N=text0/J/A=MC \"";
879 *os << "\\s(yphase_" << 1 << ") "
880 + put_bravais_type_name(latfit_tray[0].putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
881 + latfit_tray[0].putLatticeFigureOfMerit().printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
882 for (Int4 j=1, j1=2; j<isize; j++, j1++)
883 {
884 *os << "\\r\\s(yphase_" << j1 << ") "
885 + put_bravais_type_name(latfit_tray[j].putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
886 + latfit_tray[j].putLatticeFigureOfMerit().printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
887 }
888 *os << "\"\nX Legend/C/N=text0/J/A=RT/X=0.00/Y=0.00\n";
889
890 ofs.close();
891 }
Back to OSDN">Back to OSDN
 ViewVC Help
Powered by ViewVC 1.1.26  
About OSDN
Find Software
Develop Software
Help