1 /* test_core.cpp: implementation file for clipper core self-test */
2 //C Copyright (C) 2000-2006 Kevin Cowtan and University of York
3 //L
4 //L  This library is free software and is distributed under the terms
5 //L  and conditions of version 2.1 of the GNU Lesser General Public
6 //L  Licence (LGPL) with the following additional clause:
7 //L
8 //L     `You may also combine or link a "work that uses the Library" to
9 //L     produce a work containing portions of the Library, and distribute
10 //L     that work under terms of your choice, provided that you give
11 //L     prominent notice with each copy of the work that the specified
12 //L     version of the Library is used in it, and that you include or
13 //L     provide public access to the complete corresponding
14 //L     machine-readable source code for the Library including whatever
15 //L     changes were used in the work. (i.e. If you make changes to the
16 //L     Library you must distribute those, but you do not need to
17 //L     distribute source or object code to those portions of the work
18 //L     not covered by this licence.)'
19 //L
20 //L  Note that this clause grants an additional right and does not impose
21 //L  any additional restriction, and so does not affect compatibility
22 //L  with the GNU General Public Licence (GPL). If you wish to negotiate
23 //L  other terms, please contact the maintainer.
24 //L
25 //L  You can redistribute it and/or modify the library under the terms of
26 //L  the GNU Lesser General Public License as published by the Free Software
27 //L  Foundation; either version 2.1 of the License, or (at your option) any
28 //L  later version.
29 //L
30 //L  This library is distributed in the hope that it will be useful, but
31 //L  WITHOUT ANY WARRANTY; without even the implied warranty of
32 //L  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
33 //L  Lesser General Public License for more details.
34 //L
35 //L  You should have received a copy of the CCP4 licence and/or GNU
36 //L  Lesser General Public License along with this library; if not, write
37 //L  to the CCP4 Secretary, Daresbury Laboratory, Warrington WA4 4AD, UK.
38 //L  The GNU Lesser General Public can also be obtained by writing to the
39 //L  Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston,
40 //L  MA 02111-1307 USA
41 
42 
43 #include "test_core.h"
44 #include "test_data.h"
45 #include "hkl_compute.h"
46 #include "hkl_operators.h"
47 #include "resol_targetfn.h"
48 #include "atomsf.h"
49 #include "xmap.h"
50 #include "rotation.h"
51 
52 #include <algorithm>
53 
54 
55 namespace clipper {
56 
57 
operator ()()58 bool Test_core::operator() () {
59   typedef HKL_info::HKL_reference_index HRI;
60   data::Test_data data;
61   const HKL_data<data32::F_sigF>& fsig = data.hkl_data_f_sigf();
62   const HKL_data<data32::ABCD>&   abcd = data.hkl_data_abcd();
63   const Spacegroup spgr = fsig.hkl_info().spacegroup();
64   const Cell       cell = fsig.hkl_info().cell();
65   const Resolution reso = fsig.hkl_info().resolution();
66 
67   // test data sizes
68   {
69     test( "FTYPE32", sizeof(clipper::ftype32), 4 );
70     test( "ITYPE32", sizeof(clipper::itype32), 4 );
71     test( "FTYPE64", sizeof(clipper::ftype64), 8 );
72     test( "ITYPE64", sizeof(clipper::itype64), 8 );
73   }
74 
75   // test NaN functions
76   {
77     ftype64 nandd = Util::nand();
78     ftype64 nandf = Util::nanf();
79     ftype32 nanfd = Util::nand();
80     ftype32 nanff = Util::nanf();
81     test( "ISNANDD", Util::isnan( nandd )?1:0, 1 );
82     test( "ISNANDF", Util::isnan( nandf )?1:0, 1 );
83     test( "ISNANFD", Util::isnan( nanfd )?1:0, 1 );
84     test( "ISNANFF", Util::isnan( nanff )?1:0, 1 );
85     test( "IS_NANDD", Util::is_nan( nandd )?1:0, 1 );
86     test( "IS_NANDF", Util::is_nan( nandf )?1:0, 1 );
87     test( "IS_NANFD", Util::is_nan( nanfd )?1:0, 1 );
88     test( "IS_NANFF", Util::is_nan( nanff )?1:0, 1 );
89     for ( double x = -30.0; x < 30.5; x += 1.0 ) {
90       test( "ISNANX",  Util::isnan(exp(x))?1:0, 0 );
91       test( "IS_NANX",  Util::is_nan(exp(x))?1:0, 0 );
92     }
93   }
94 
95   // test Sim functions for consistency
96   for ( double x=0.1; x<50.0; x*=1.5 ) {
97     test( "SIMINV", x, Util::invsim(Util::sim(x)), 0.001*x );
98     test( "SIMNEG", Util::sim(x), -Util::sim(-x), 0.001 );
99     test( "SIMINT", Util::sim(x), (Util::sim_integ(x+0.0001)-Util::sim_integ(x))/0.0001, 0.001 );
100     test( "SIMDRV", Util::sim_deriv(x), (Util::sim(x+0.0001)-Util::sim(x))/0.0001, 0.001 );
101     if ( x < 2.0 ) test( "SIM_I0", Util::sim(x), (Util::bessel_i0(x+0.0001)/Util::bessel_i0(x)-1.0)/0.0001, 0.001 );
102   }
103 
104   // test HL<->abcd conversion
105   for ( HRI ih = abcd.first(); !ih.last(); ih.next() ) {
106     data32::Compute_phifom_from_abcd cap;
107     data32::Compute_abcd_from_phifom cpa;
108     data32::Phi_fom phiw1, phiw2;
109     data32::ABCD    abcd1, abcd2;
110     LogPhaseProb<24> q(ih.hkl_class());
111     q.set_abcd( abcd[ih] ); // 1st route
112     q.get_phi_fom( phiw1 );
113     abcd1 = cpa( ih, phiw1 );
114     phiw1 = cap( ih, abcd1 );
115     phiw2 = cap( ih, abcd[ih] ); // 2nd route
116     q.set_phi_fom( phiw2 );
117     q.get_abcd( abcd2 );
118     phiw2 = cap( ih, abcd2 );
119     test( "HLW", phiw1.fom(), phiw2.fom(), 0.003 );
120   }
121 
122   // test map calculation
123   {
124     HKL_data<data32::Phi_fom> pw( fsig.hkl_info() );
125     HKL_data<data32::F_phi> fp1( fsig.hkl_info() ), fp2( fsig.hkl_info() );
126     pw.compute( abcd, data32::Compute_phifom_from_abcd() );
127     fp1.compute( fsig, pw, data32::Compute_fphi_from_fsigf_phifom() );
128     Grid_sampling grid( spgr, cell, reso );
129     Xmap<float> xmap( spgr, cell, grid );
130     xmap.fft_from( fp1 );
131     xmap.fft_to( fp2 );
132     for ( HRI ih = fp1.first(); !ih.last(); ih.next() ) {
133       std::complex<float> ab1(0.0,0.0), ab2(0.0,0.0);
134       if ( !fp1[ih].missing() ) ab1 = fp1[ih];
135       if ( !fp2[ih].missing() ) ab2 = fp2[ih];
136       test( "FFT-A", ab1.real(), ab2.real(), 0.01 );
137       test( "FFT-B", ab1.imag(), ab2.imag(), 0.01 );
138     }
139   }
140 
141   // test resolution ordinals
142   {
143     std::vector<clipper::ftype> resols;
144     for ( HRI ih = fsig.first(); !ih.last(); ih.next() )
145       resols.push_back( ih.invresolsq() );
146     clipper::Generic_ordinal ordinal;
147     ordinal.init(resols);
148     clipper::Generic_ordinal ordinv = ordinal;
149     ordinv.invert();
150     clipper::Resolution_ordinal resord;
151     resord.init( fsig.hkl_info(), 1.0 );
152     std::sort( resols.begin(), resols.end() );
153     for ( int i = 0; i < resols.size(); i++ )
154       test( "RESORD", resols[i], ordinv.ordinal(resord.ordinal(resols[i])), 0.001 );
155   }
156 
157   // test resolution functions
158   {
159     std::vector<ftype> param( 10, 1.0 );
160     BasisFn_spline basisfn( fsig, 10 );
161     TargetFn_meanFnth<data32::F_sigF> targetfn( fsig, 2.0 );
162     ResolutionFn rfn( fsig.hkl_info(), basisfn, targetfn, param );
163     test( "RESFN0", rfn.params()[0], 229690.9746, 0.1 );
164     test( "RESFN1", rfn.params()[1], 216481.7609, 0.1 );
165     test( "RESFN2", rfn.params()[2],  78484.9498, 0.1 );
166     test( "RESFN3", rfn.params()[3], 148774.2654, 0.1 );
167     test( "RESFN4", rfn.params()[4],  69255.6000, 0.1 );
168     test( "RESFN5", rfn.params()[5], 143032.5088, 0.1 );
169     test( "RESFN6", rfn.params()[6], 110371.3524, 0.1 );
170     test( "RESFN7", rfn.params()[7], 108711.3487, 0.1 );
171     test( "RESFN8", rfn.params()[8], 150487.5496, 0.1 );
172     test( "RESFN9", rfn.params()[9], 141713.7420, 0.1 );
173   }
174 
175   // test atom shape function derivatives
176   {
177     ftype d = 0.001;
178     Coord_orth co( 1.0, 2.0, 3.0 );
179     AtomShapeFn sf( co, "N", 0.25, 1.0 );
180     AtomShapeFn sfx( co+Coord_orth(d,0.0,0.0), "N", 0.25, 1.0 );
181     AtomShapeFn sfy( co+Coord_orth(0.0,d,0.0), "N", 0.25, 1.0 );
182     AtomShapeFn sfz( co+Coord_orth(0.0,0.0,d), "N", 0.25, 1.0 );
183     AtomShapeFn sfo( co, "N", 0.25, 1.0+d );
184     AtomShapeFn sfu( co, "N", 0.251, 1.0 );
185     AtomShapeFn sfx2( co+Coord_orth(2.0*d,0.0,0.0), "N", 0.25, 1.0 );
186     AtomShapeFn sfy2( co+Coord_orth(0.0,2.0*d,0.0), "N", 0.25, 1.0 );
187     AtomShapeFn sfz2( co+Coord_orth(0.0,0.0,2.0*d), "N", 0.25, 1.0 );
188     AtomShapeFn sfo2( co, "N", 0.25, 1.0+d+d );
189     AtomShapeFn sfu2( co, "N", 0.252, 1.0 );
190     AtomShapeFn sfxy( co+Coord_orth(d,d,0.0), "N", 0.25, 1.0 );
191     AtomShapeFn sfyz( co+Coord_orth(0.0,d,d), "N", 0.25, 1.0 );
192     AtomShapeFn sfzx( co+Coord_orth(d,0.0,d), "N", 0.25, 1.0 );
193     U_aniso_orth uai( 0.25, 0.25, 0.25, 0.0, 0.0, 0.0 );
194     AtomShapeFn sfuai( co, "N", uai, 1.0 );
195     std::vector<AtomShapeFn::TYPE> params;
196     params.push_back( AtomShapeFn::X );
197     params.push_back( AtomShapeFn::Y );
198     params.push_back( AtomShapeFn::Z );
199     params.push_back( AtomShapeFn::Occ );
200     params.push_back( AtomShapeFn::Uiso );
201     params.push_back( AtomShapeFn::U11 );
202     params.push_back( AtomShapeFn::U22 );
203     params.push_back( AtomShapeFn::U33 );
204     params.push_back( AtomShapeFn::U12 );
205     params.push_back( AtomShapeFn::U13 );
206     params.push_back( AtomShapeFn::U23 );
207     sf.agarwal_params() = params;
208     for ( int i = 0; i < 100; i++ ) {
209       Coord_orth c2 = co + Coord_orth( 0.11*(i%5-1.9), 0.13*(i%7-2.8), 0.15*(i%9-3.7) );
210       double f;
211       std::vector<double> g(11);
212       Matrix<double> c(11,11);
213       test( "ATOMSF-A", sfuai.rho(c2), sf.rho(c2), 1.0e-8 );
214       sf.rho_curv( c2, f, g, c );
215       test( "ATOMSF-G", (sfx.rho(c2)-sf.rho(c2))/d, g[0], 0.01 );
216       test( "ATOMSF-G", (sfy.rho(c2)-sf.rho(c2))/d, g[1], 0.01 );
217       test( "ATOMSF-G", (sfz.rho(c2)-sf.rho(c2))/d, g[2], 0.01 );
218       test( "ATOMSF-G", (sfo.rho(c2)-sf.rho(c2))/d, g[3], 0.01 );
219       test( "ATOMSF-G", (sfu.rho(c2)-sf.rho(c2))/d, g[4], 0.05 );
220       test( "ATOMSF-C", (sfx2.rho(c2)-2*sfx.rho(c2)+sf.rho(c2))/(d*d), c(0,0), 0.1 );
221       test( "ATOMSF-C", (sfy2.rho(c2)-2*sfy.rho(c2)+sf.rho(c2))/(d*d), c(1,1), 0.1 );
222       test( "ATOMSF-C", (sfz2.rho(c2)-2*sfz.rho(c2)+sf.rho(c2))/(d*d), c(2,2), 0.1 );
223       test( "ATOMSF-C", (sfxy.rho(c2)-sfx.rho(c2)-sfy.rho(c2)+sf.rho(c2))/(d*d), c(0,1), 0.1 );
224       test( "ATOMSF-C", (sfyz.rho(c2)-sfy.rho(c2)-sfz.rho(c2)+sf.rho(c2))/(d*d), c(1,2), 0.1 );
225       test( "ATOMSF-C", (sfzx.rho(c2)-sfz.rho(c2)-sfx.rho(c2)+sf.rho(c2))/(d*d), c(2,0), 0.1 );
226     }
227     for ( int j = 0; j < 20; j++ ) {
228       ftype x = 0.19*(j%5-1.9);
229       ftype y = 0.15*(j%7-2.8);
230       ftype z = 0.13*(j%9-3.7);
231       U_aniso_orth ua  ( x*x+0.2, y*y+0.2, z*z+0.2, y*z, z*x, x*y );
232       U_aniso_orth ua00( ua.mat00()+d, ua.mat11(), ua.mat22(),
233 			 ua.mat01(), ua.mat02(), ua.mat12() );
234       U_aniso_orth ua11( ua.mat00(), ua.mat11()+d, ua.mat22(),
235 			 ua.mat01(), ua.mat02(), ua.mat12() );
236       U_aniso_orth ua22( ua.mat00(), ua.mat11(), ua.mat22()+d,
237 			 ua.mat01(), ua.mat02(), ua.mat12() );
238       U_aniso_orth ua01( ua.mat00(), ua.mat11(), ua.mat22(),
239 			 ua.mat01()+d, ua.mat02(), ua.mat12() );
240       U_aniso_orth ua02( ua.mat00(), ua.mat11(), ua.mat22(),
241 			 ua.mat01(), ua.mat02()+d, ua.mat12() );
242       U_aniso_orth ua12( ua.mat00(), ua.mat11(), ua.mat22(),
243 			 ua.mat01(), ua.mat02(), ua.mat12()+d );
244       AtomShapeFn sfua ( co, "N", ua, 1.0 );
245       AtomShapeFn sfuax( co+Coord_orth(d,0.0,0.0), "N", ua, 1.0 );
246       AtomShapeFn sfuay( co+Coord_orth(0.0,d,0.0), "N", ua, 1.0 );
247       AtomShapeFn sfuaz( co+Coord_orth(0.0,0.0,d), "N", ua, 1.0 );
248       AtomShapeFn sfuao( co, "N", ua, 1.0+d );
249       AtomShapeFn sfua00( co, "N", ua00, 1.0 );
250       AtomShapeFn sfua11( co, "N", ua11, 1.0 );
251       AtomShapeFn sfua22( co, "N", ua22, 1.0 );
252       AtomShapeFn sfua01( co, "N", ua01, 1.0 );
253       AtomShapeFn sfua02( co, "N", ua02, 1.0 );
254       AtomShapeFn sfua12( co, "N", ua12, 1.0 );
255       sfua.agarwal_params() = params;
256       for ( int i = 0; i < 50; i++ ) {
257 	Coord_orth c2 = co + Coord_orth( 0.11*(i%5-1.9), 0.13*(i%7-2.8), 0.15*(i%9-3.7) );
258 	double f;
259 	std::vector<double> g(11);
260 	Matrix<double> c(11,11);
261 	sfua.rho_grad( c2, f, g );
262 	test( "ATOMSF-AG", ( sfuax.rho(c2)-sfua.rho(c2))/d,  g[0], 0.01 );
263 	test( "ATOMSF-AG", ( sfuay.rho(c2)-sfua.rho(c2))/d,  g[1], 0.01 );
264 	test( "ATOMSF-AG", ( sfuaz.rho(c2)-sfua.rho(c2))/d,  g[2], 0.01 );
265 	test( "ATOMSF-AG", ( sfuao.rho(c2)-sfua.rho(c2))/d,  g[3], 0.01 );
266 	test( "ATOMSF-AG", (sfua00.rho(c2)-sfua.rho(c2))/d,  g[5], 0.05 );
267 	test( "ATOMSF-AG", (sfua11.rho(c2)-sfua.rho(c2))/d,  g[6], 0.05 );
268 	test( "ATOMSF-AG", (sfua22.rho(c2)-sfua.rho(c2))/d,  g[7], 0.05 );
269 	test( "ATOMSF-AG", (sfua01.rho(c2)-sfua.rho(c2))/d,  g[8], 0.05 );
270 	test( "ATOMSF-AG", (sfua02.rho(c2)-sfua.rho(c2))/d,  g[9], 0.05 );
271 	test( "ATOMSF-AG", (sfua12.rho(c2)-sfua.rho(c2))/d, g[10], 0.05 );
272       }
273     }
274   }
275 
276   // test spacegroups
277   {
278     const char *pgs[] = {"-P 1", "-P 2", "-P 2y", "-P 2x", "-P 2\"", "-P 2y\"", "-P 2x\"", "-P 2'", "-P2y'", "-P 2x'", "-P 2 2", "-P 2 2\"", "-P 2 2\"(y,z,x)", "-P 2 2\"(z,x,y)", "-P3", "-P 3 (y,z,x)", "-P 3 (z,x,y)", "-P 3 (-x,y,z)", "-P 3 (y,z,-x)", "-P 3 (z,-x,y)", "-P 3*", "-P 3* (-x,y,z)", "-P 3* (x,-y,z)", "-P 3* (x,y,-z)", "-P 3 2", "-P 3 2 (y,z,x)", "-P 3 2 (z,x,y)", "-P 3* 2", "-P 3* 2 (-x,y,z)", "-P 3* 2 (x,-y,z)", "-P 3* 2 (-x,-y,z)", "-P 3 2\"", "-P 3 2\"(z,x,y)", "-P 3 2\"(y,z,x)", "-P 3 2\"(-x,y,z)", "-P 3 2\"(z,-x,y)", "-P 3 2\"(y,z,-x)", "-P 4", "-P 4 (y,z,x)", "-P 4 (z,x,y)", "-P 4 2", "-P 4 2 (y,z,x)", "-P 4 2 (z,x,y)", "-P 6", "-P 6 (y,z,x)", "-P 6 (z,x,y)", "-P 6 2", "-P 6 2 (y,z,x)", "-P 6 2 (z,x,y)", "-P 2 2 3", "-P 4 2 3" };
279     std::vector<String> hallsymbols;
280     for ( int i = 0; i < data::sgdata_size; i++ )
281       hallsymbols.push_back( data::sgdata[i].hall );
282     for ( int i = 0; i < sizeof(pgs)/sizeof(pgs[0]); i++ )
283       hallsymbols.push_back( pgs[i] );
284     Cell cellc( Cell_descr( 37, 37, 37 ) );
285     Cell cellha( Cell_descr( 37, 37, 37, 120, 90, 90 ) );
286     Cell cellhb( Cell_descr( 37, 37, 37, 90, 120, 90 ) );
287     Cell cellhc( Cell_descr( 37, 37, 37, 90, 90, 120 ) );
288     Cell cellha1( Cell_descr( 37, 37, 37, 60, 90, 90 ) );
289     Cell cellhb1( Cell_descr( 37, 37, 37, 90, 60, 90 ) );
290     Cell cellhc1( Cell_descr( 37, 37, 37, 90, 90, 60 ) );
291     Cell cell;
292     String symbol;
293     Spacegroup sg;
294     for ( int s = 0; s < hallsymbols.size(); s++ ) {
295       try {
296 	symbol = hallsymbols[s];
297 	sg = Spacegroup( Spgr_descr( symbol, Spgr_descr::Hall ) );
298 	// identify trigonal/hexagonal groups
299 	cell = cellc;
300 	for ( int sym = 1; sym < sg.num_symops(); sym++ ) {
301 	  if ( ( sg.symop(sym).rot()(1,1) * sg.symop(sym).rot()(1,2) == -1 ) ||
302 	       ( sg.symop(sym).rot()(2,1) * sg.symop(sym).rot()(2,2) == -1 ) )
303 	    cell = cellha;
304 	  if ( ( sg.symop(sym).rot()(0,0) * sg.symop(sym).rot()(0,2) == -1 ) ||
305 	       ( sg.symop(sym).rot()(2,0) * sg.symop(sym).rot()(2,2) == -1 ) )
306 	    cell = cellhb;
307 	  if ( ( sg.symop(sym).rot()(0,0) * sg.symop(sym).rot()(0,1) == -1 ) ||
308 	       ( sg.symop(sym).rot()(1,0) * sg.symop(sym).rot()(1,1) == -1 ) )
309 	    cell = cellhc;
310 	  if ( ( sg.symop(sym).rot()(1,1) * sg.symop(sym).rot()(1,2) == 1 ) ||
311 	       ( sg.symop(sym).rot()(2,1) * sg.symop(sym).rot()(2,2) == 1 ) )
312 	    cell = cellha1;
313 	  if ( ( sg.symop(sym).rot()(0,0) * sg.symop(sym).rot()(0,2) == 1 ) ||
314 	       ( sg.symop(sym).rot()(2,0) * sg.symop(sym).rot()(2,2) == 1 ) )
315 	    cell = cellhb1;
316 	  if ( ( sg.symop(sym).rot()(0,0) * sg.symop(sym).rot()(0,1) == 1 ) ||
317 	       ( sg.symop(sym).rot()(1,0) * sg.symop(sym).rot()(1,1) == 1 ) )
318 	    cell = cellhc1;
319 	}
320 	for ( int i = 0; i < 100; i++ ) {
321 	  HKL rfl( i%5-2, i%7-3, i%9-4 );
322 	  ftype s0 = rfl.invresolsq(cell);
323 	  for ( int sym = 1; sym < sg.num_symops(); sym++ ) {
324 	    ftype s1 = rfl.transform(sg.symop(sym)).invresolsq(cell);
325 	    test( "SG "+symbol, s0, s1, 1.0e-12 );
326 	  }
327 	}
328 	Grid_sampling grid( sg, cell, Resolution( 4.0 ) );
329 	Xmap<ftype32> xmap( sg, cell, grid );
330       } catch ( Message_base ) {
331 	test( "SG "+symbol, sg.spacegroup_number(), -1 );
332       }
333     }
334   }
335 
336   // test rotations
337   {
338     for ( ftype x = -1.0; x < 1.0; x += 0.02 )
339       for ( ftype y = -1.0; y < 1.0; y += 0.02 )
340 	for ( ftype z = -1.0; z < 1.0; z += 0.02 ) {
341 	  ftype s = x*x + y*y + z*z;
342 	  if ( s < 1.0 ) {
343 	    ftype w = sqrt( 1.0 - s );
344 	    Rotation rot( w, x, y, z );
345 	    Rotation rotinv = rot.inverse();
346 	    Rotation r1( rot.matrix() );
347 	    Rotation r2( rot.polar_ccp4() );
348 	    test( "ROT/MAT "+rot.format(), (rotinv*r1).abs_angle(), 0.0, 1.0e-6 );
349 	    test( "ROT/POL "+rot.format(), (rotinv*r2).abs_angle(), 0.0, 1.0e-6 );
350 	  }
351 	}
352     /*
353     Mat33<> mat1( -0.18332,  0.02511, -0.98273,
354 		  0.02184, -0.99932, -0.02960,
355 		 -0.98281, -0.02689,  0.18265 );
356     Rotation r1( mat1 );
357     Mat33<> mat2( r1.matrix() );
358     Rotation r2( mat2 );
359     std::cout << mat1.format() << std::endl << r1.format() << std::endl << mat2.format() << std::endl << r2.format() << std::endl << std::endl;
360     Rotation r1a( mat1 ); r1a.norm();
361     Mat33<> mat2a( r1a.matrix() );
362     Rotation r2a( mat2a );
363     Mat33<> mat3a( r2a.matrix() );
364     Rotation r3a( mat3a );
365     std::cout << mat1.format() << mat1.det() << std::endl << r1a.format() << std::endl << mat2a.format() << mat2a.det() << std::endl << r2a.format() << std::endl << mat3a.format() << mat3a.det() << std::endl << r3a.format() << std::endl << std::endl;
366     Rotation r3(0.02,0.42,0.02,-0.50);
367     std::cout << r3.format() << std::endl << Rotation(r3.matrix()).format() << std::endl;
368     */
369   }
370 
371   return ( error_count == 0 );
372 }
373 
374 
375 } // namespace clipper
376