Develop and Download Open Source Software

Powder auto-indexing software Conograph

Browse Subversion Repository

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

Annotation of /Conograph/trunk/src/p_out_indexing.cc

Parent Directory Parent Directory | Revision Log Revision Log

Revision 25 - (hide 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 rtomiyasu 3 /*
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 rtomiyasu 25 #include "lattice_symmetry/VCLatticeFigureOfMeritToCheckSymmetry.hh"
32 rtomiyasu 3 #include "lattice_symmetry/LatticeFigureOfMeritToDisplay.hh"
33     #include "utility_func/covar_matrix.hh"
34     #include "utility_func/zmath.hh"
35 rtomiyasu 25 #include "utility_lattice_reduction/matrix_NbyN.hh"
36 rtomiyasu 3 #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 rtomiyasu 25 #include "utility_data_structure/VecDat3.hh"
41 rtomiyasu 3 #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 rtomiyasu 25 void printHKLdata(const vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[],
331 rtomiyasu 3 const OutputInfo outinfo[],
332     const eSortCriterion& sort_criterion,
333     const ControlParam& cdata,
334     const PeakPosData& pdata,
335     const string& fname)
336     {
337 rtomiyasu 25 static const Int4 NUM_LS = put_number_of_bravais_types();
338 rtomiyasu 3
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 rtomiyasu 25 for(Int4 i=0; i<putNumberOfSortCriterion(); i++)
385 rtomiyasu 3 {
386     *os << ", " << putLabel(eSortCriterion(i));
387     }
388     *os << ", VOL\n";
389    
390 rtomiyasu 25 for(Int4 k=NUM_LS-1; k>=0; k--)
391 rtomiyasu 3 {
392     os->width(17);
393 rtomiyasu 25 *os << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << " : ";
394 rtomiyasu 3 os->width(7);
395     *os << outinfo[k].putCandidateNumToOutput();
396    
397 rtomiyasu 25 const Int4 index = outinfo[k].putIndex( outinfo[k].putLatticeLabelSelectedFromListToOutput(sort_criterion) );
398 rtomiyasu 3 if( index >= 0 )
399     {
400     const LatticeFigureOfMerit& lat_fom = lattice_result[k][index].putLatticeFigureOfMerit();
401 rtomiyasu 25 const LatticeFigureOfMerit::SetOfFigureOfMerit& set_fom = lat_fom.putFiguresOfMerit();
402 rtomiyasu 3
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 rtomiyasu 25 for(Int4 k=NUM_LS-1; k>=0; k--)
437 rtomiyasu 3 {
438     os->width(18);
439 rtomiyasu 25 *os << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << ", Best Score : Lattice constants, label.\n";
440 rtomiyasu 3 // if( lattice_result[k].empty() )
441     // {
442     // *os << "\n";
443     // continue;
444     // }
445 rtomiyasu 25 for(Int4 i=0; i<putNumberOfSortCriterion(); i++)
446 rtomiyasu 3 {
447     os->width(15);
448     *os << putLabel(eSortCriterion(i)) << " = ";
449    
450 rtomiyasu 25 const Int4 index = outinfo[k].putIndex( outinfo[k].putLatticeLabelSelectedAmongAll(eSortCriterion(i)) );
451 rtomiyasu 3 os->width(12);
452     if( index >= 0 )
453     {
454     const LatticeFigureOfMerit& lat_fom = lattice_result[k][index].putLatticeFigureOfMerit();
455 rtomiyasu 25 const LatticeFigureOfMerit::SetOfFigureOfMerit& set_fom = lat_fom.putFiguresOfMerit();
456 rtomiyasu 3
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 rtomiyasu 25 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 rtomiyasu 3
521     if( selected_lattice_index >= 0 )
522     {
523 rtomiyasu 25 const VCLatticeFigureOfMeritToCheckSymmetry& selected_lattice = lattice_result[selected_lattice_info.first][selected_lattice_index];
524 rtomiyasu 3
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 rtomiyasu 25 selected_lattice.printLatticeInformation(lattice_result, outinfo,
533     cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
534     cdata.putResolution(), label_start, os);
535 rtomiyasu 3
536     label_start--;
537     os->width(label_start);
538     *os << "" << "</SelectedLatticeCandidate>\n\n";
539     }
540    
541 rtomiyasu 25 for(Int4 k=NUM_LS-1; k>=0; k--)
542 rtomiyasu 3 {
543     *os << "\n";
544     os->width(label_start);
545 rtomiyasu 25 *os << "" << "<!-- Candidates for " << put_bravais_type_name(eBravaisType(k), cdata.putBaseCenteredAxis()) << " -->\n\n";
546 rtomiyasu 3
547     const Int4 num_topo = lattice_result[k].size();
548    
549     for(Int4 n=0; n<num_topo; n++)
550     {
551 rtomiyasu 25 const VCLatticeFigureOfMeritToCheckSymmetry& ans = lattice_result[k][n];
552 rtomiyasu 3 if( !outinfo[k].IsOutput(ans.putLatticeLabel()) ) continue;
553    
554     os->width(label_start);
555     *os << "" << "<LatticeCandidate number=\"" << ans.putStringLabel() << "\">\n";
556     label_start++;
557 rtomiyasu 25 ans.printLatticeInformation(lattice_result, outinfo, cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
558     0.06, label_start, os);
559 rtomiyasu 3
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 rtomiyasu 25 ans.printLatticeInformation(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(),
615     cdata.putResolution(), label_start, os);
616 rtomiyasu 3 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 rtomiyasu 25
636 rtomiyasu 3 void printPeakPosition(
637     const ControlParam& cdata,
638     const PeakPosData& pdata,
639     const LatticeFigureOfMeritToDisplay& latfit,
640 rtomiyasu 25 const vector<LatticeFigureOfMerit>& lattices_same_q,
641 rtomiyasu 3 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 rtomiyasu 25 *os << "WAVES/O dphase_0, xphase_0, yphase_0, h_0, k_0, l_0\n";
652 rtomiyasu 3
653 rtomiyasu 25 *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 rtomiyasu 3 {
698 rtomiyasu 25 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 rtomiyasu 3 *os << endl;
710    
711     *os << "BEGIN\n";
712    
713 rtomiyasu 25 const vector<HKL_Q>& cal_hkl_tray = latfit2.putCalMillerIndices();
714 rtomiyasu 3 Vec_DP cal_pos_tray;
715 rtomiyasu 25 latfit2.putCalculatedPeakPosInRange(cdata, cal_pos_tray);
716 rtomiyasu 3
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 rtomiyasu 25 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 rtomiyasu 3 {
767 rtomiyasu 25 *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 rtomiyasu 3 }
772    
773     *os << "X Legend/C/N=text0/J/A=MC \"";
774 rtomiyasu 25 *os << "\\s(yphase_" << 0 << ") "
775     + put_bravais_type_name(latfit.putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
776 rtomiyasu 3 + latfit.putLatticeFigureOfMerit().printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
777 rtomiyasu 25 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 rtomiyasu 3
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 rtomiyasu 25 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 rtomiyasu 3 {
873     *os << "X AppendToGraph yphase_" << j1 << " vs xphase_" << j1 << endl;
874 rtomiyasu 16 *os << "X ModifyGraph offset(yphase_" << j1 << ")={0," << offset << "},mode(yphase_" << j1 << ")=3,marker(yphase_" << j1;
875 rtomiyasu 3 *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 rtomiyasu 25 + put_bravais_type_name(latfit_tray[0].putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
881 rtomiyasu 3 + 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 rtomiyasu 25 + put_bravais_type_name(latfit_tray[j].putLatticeFigureOfMerit().enumBravaisType(), cdata.putBaseCenteredAxis()) + " "
886 rtomiyasu 3 + latfit_tray[j].putLatticeFigureOfMerit().printOptimizedLatticeConstants(cdata.putBaseCenteredAxis(), cdata.putRhombohedralAxis(), 3);
887     }
888 rtomiyasu 25 *os << "\"\nX Legend/C/N=text0/J/A=RT/X=0.00/Y=0.00\n";
889 rtomiyasu 3
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