EMAN2
symmetry.cpp
Go to the documentation of this file.
1/*
2 * Author: David Woolford, 09/23/2008 (woolford@bcm.edu)
3 * Copyright (c) 2000-2006 Baylor College of Medicine
4 *
5 * This software is issued under a joint BSD/GNU license. You may use the
6 * source code in this file under either license. However, note that the
7 * complete EMAN2 and SPARX software packages have some GPL dependencies,
8 * so you are responsible for compliance with the licenses of these packages
9 * if you opt to use BSD licensing. The warranty disclaimer below holds
10 * in either instance.
11 *
12 * This complete copyright notice must be included in any revised version of the
13 * source code. Additional authorship citations may be added, but existingx
14 * author citations must be preserved.
15 *
16 * This program is free software; you can redistribute it and/or modify
17 * it under the terms of the GNU General Public License as published by
18 * the Free Software Foundation; either version 2 of the License, or
19 * (at your option) any later version.
20 *
21 * This program is distributed in the hope that it will be useful,
22 * but WITHOUT ANY WARRANTY; without even the implied warranty of
23 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 * GNU General Public License for more details.
25 *
26 * You should have received a copy of the GNU General Public License
27 * along with this program; if not, write to the Free Software
28 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
29 *
30 * */
31
32#include "symmetry.h"
33#include "transform.h"
34#include "vec3.h"
35#include "exception.h"
36#include "util.h"
37
38using namespace EMAN;
39
40const string CSym::NAME = "c";
41const string DSym::NAME = "d";
42const string HSym::NAME = "h";
43const string TetrahedralSym::NAME = "tet";
44const string OctahedralSym::NAME = "oct";
45const string IcosahedralSym::NAME = "icos";
46const string Icosahedral2Sym::NAME = "icos2";
47const string EmanOrientationGenerator::NAME = "eman";
48const string SingleOrientationGenerator::NAME = "single";
49const string SaffOrientationGenerator::NAME = "saff";
50const string EvenOrientationGenerator::NAME = "even";
51const string RandomOrientationGenerator::NAME = "rand";
52const string OptimumOrientationGenerator::NAME = "opt";
53
54template <> Factory < Symmetry3D >::Factory()
55{
56 force_add<CSym>();
57 force_add<DSym>();
58 force_add<HSym>();
59 force_add<TetrahedralSym>();
60 force_add<OctahedralSym>();
61 force_add<IcosahedralSym>();
62 force_add<Icosahedral2Sym>();
63}
64
65void EMAN::dump_symmetries()
66{
67 dump_factory < Symmetry3D > ();
68}
69
70map<string, vector<string> > EMAN::dump_symmetries_list()
71{
72 return dump_factory_list < Symmetry3D > ();
73}
74
75template <>
76Symmetry3D* Factory < Symmetry3D >::get(const string & instancename_)
77{
78 init();
79
80 string instancename = Util::str_to_lower(instancename_);
81
82 unsigned int n = instancename.size();
83 if ( n == 0 ) throw NotExistingObjectException(instancename, "Empty instance name!");
84
85 char leadingchar = instancename[0];
86 if (leadingchar == 'c' || leadingchar == 'd' || leadingchar == 'h' ) {
87 Dict parms;
88 if (n > 1) {
89 int nsym = atoi(instancename.c_str() + 1);
90 parms["nsym"] = nsym;
91 }
92
93 if (leadingchar == 'c') {
94 return get("c",parms);
95 }
96 if (leadingchar == 'd') {
97 return get("d",parms);
98 }
99 if (leadingchar == 'h') {
100 int nstart=1,nsym=1;
101 float daz=0.,tz=0.,maxtilt=90.0;
102 string temp;
103 temp=instancename;
104 temp.erase(0,1);
105 if (sscanf(temp.c_str(),"%d:%d:%f:%f:%f",&nsym,&nstart,&daz,&tz,&maxtilt)<4) {
106 sscanf(temp.c_str(),"%d,%d,%f,%f,%f",&nsym,&nstart,&daz,&tz,&maxtilt);
107 }
108 parms["nstart"]=nstart;
109 parms["nsym"]=nsym;
110 parms["daz"]=daz;
111 parms["tz"]=tz;
112 parms["maxtilt"]=maxtilt;
113 return get("h",parms);
114 }
115
116// delete lc;
117 }
118 else if ( instancename == "icos" || instancename == "oct" || instancename == "tet" || instancename == "icos2" || instancename == "icos5" )
119 {
120 if (instancename == "icos5") {
121 instancename = "icos";
122 }
123
124 map < string, InstanceType >::iterator fi =
125 my_instance->my_dict.find(instancename);
126 if (fi != my_instance->my_dict.end()) {
127 return my_instance->my_dict[instancename] ();
128 }
129
130 string lower = instancename;
131 for (unsigned int i=0; i<lower.length(); i++) lower[i]=tolower(lower[i]);
132
133 fi = my_instance->my_dict.find(lower);
134 if (fi != my_instance->my_dict.end()) {
135 return my_instance->my_dict[lower] ();
136 }
137
138 throw NotExistingObjectException(instancename, "No such an instance existing");
139 }
140 else throw NotExistingObjectException(instancename, "No such an instance existing");
141
142 throw NotExistingObjectException(instancename, "No such an instance existing");
143}
144
145template <> Factory < OrientationGenerator >::Factory()
146{
147 force_add<EmanOrientationGenerator>();
148 force_add<SingleOrientationGenerator>();
149 force_add<RandomOrientationGenerator>();
150 force_add<EvenOrientationGenerator>();
151 force_add<SaffOrientationGenerator>();
152 force_add<OptimumOrientationGenerator>();
153}
154
155
156
157void EMAN::dump_orientgens()
158{
159 dump_factory < OrientationGenerator > ();
160}
161
162map<string, vector<string> > EMAN::dump_orientgens_list()
163{
164 return dump_factory_list < OrientationGenerator > ();
165}
166
167vector<Transform> Symmetry3D::gen_orientations(const string& generatorname, const Dict& parms)
168{
169 ENTERFUNC;
170 vector<Transform> ret;
171 OrientationGenerator *g = Factory < OrientationGenerator >::get(Util::str_to_lower(generatorname), parms);
172 if (g) {
173 ret = g->gen_orientations(this);
174 if( g )
175 {
176 delete g;
177 g = 0;
178 }
179 }
180 else throw;
181
182 EXITFUNC;
183
184 return ret;
185}
186
187void OrientationGenerator::get_az_max(const Symmetry3D* const sym, const float& altmax, const bool inc_mirror, const float& alt_iterator,const float& h,bool& d_odd_mirror_flag, float& azmax_adjusted) const
188{
189
190 if ( sym->is_d_sym() && alt_iterator == altmax && ( (sym->get_nsym())/2 % 2 == 1 )) {
191 if (inc_mirror) {
192 azmax_adjusted /= 4.0f;
193 d_odd_mirror_flag = true;
194 }
195 else azmax_adjusted /= 2.0f;
196 }
197 else if (sym->is_d_sym() && alt_iterator == altmax && ( (sym->get_nsym())/2 % 2 == 0 ) && inc_mirror) {
198 azmax_adjusted /= 2.0f;
199 }
200 // if this is odd c symmetry, and we're at the equator, and we're excluding the mirror then
201 // half the equator is redundant (it is the mirror of the other half)
202 else if (sym->is_c_sym() && !inc_mirror && alt_iterator == altmax && (sym->get_nsym() % 2 == 1 ) ){
203 azmax_adjusted /= 2.0f;
204 }
205 // at the azimuthal boundary in c symmetry and tetrahedral symmetry we have come
206 // full circle, we must not include it
207 else if (sym->is_c_sym() || sym->is_tet_sym() ) {
208 azmax_adjusted -= h/4.0f;
209 }
210 // If we're including the mirror then in d and icos and oct symmetry the azimuthal
211 // boundary represents coming full circle, so must be careful to exclude it
212 else if (inc_mirror && ( sym->is_d_sym() || sym->is_platonic_sym() ) ) {
213 azmax_adjusted -= h/4.0f;
214 }
215 // else do nothing - this means that we're including the great arc traversing
216 // the full range of permissable altitude angles at azmax.
217 // This happens in d symmetry, and in the icos and oct symmetries, when the mirror
218 // portion of the asymmetric unit is being excluded
219
220}
221
222
223//vector<Vec3f> OrientationGenerator::get_graph_opt(const Symmetry3D* const sym) const {
224// bool inc_mirror = params.set_default("inc_mirror",false);
225// vector<Vec3f> seeds = get_asym_unit_points(inc_mirror);
226// vector<Vec3i> connections;
227// if (seeds.size() == 3) {
228// connections.push_back(Vec2i(0,1,2));
229// } else {
230// throw;
231//
232// }
233//}
234
235
236float OrientationGenerator::get_optimal_delta(const Symmetry3D* const sym, const int& n) const
237{
238
239// float delta_soln = 360.0f/sym->get_max_csym();
240 float delta_soln = 180.0f;
241 float delta_upper_bound = delta_soln;
242 float delta_lower_bound = 0.0f;
243
244 int prev_tally = -1;
245 // This is an example of a divide and conquer approach, the possible values of delta are searched
246 // like a binary tree
247
248 bool soln_found = false;
249
250 while ( soln_found == false ) {
251 int tally = get_orientations_tally(sym,delta_soln);
252 if ( tally == n ) soln_found = true;
253 else if ( (delta_upper_bound - delta_lower_bound) < 0.0001 ) {
254 // If this is the case, the requested number of projections is practically infeasible
255 // in which case just return the nearest guess
256 soln_found = true;
257 delta_soln = (delta_upper_bound+delta_lower_bound)/2.0f;
258 }
259 else if (tally < n) {
260 delta_upper_bound = delta_soln;
261 delta_soln = delta_soln - (delta_soln-delta_lower_bound)/2.0f;
262 }
263 else /* tally > n*/{
264 delta_lower_bound = delta_soln;
265 delta_soln = delta_soln + (delta_upper_bound-delta_soln)/2.0f;
266 }
267 prev_tally = tally;
268 }
269
270 return delta_soln;
271}
272
273bool OrientationGenerator::add_orientation(vector<Transform>& v, const float& az, const float& alt) const
274{
275 bool randphi = params.set_default("random_phi",false);
276 float phi = 0.0f;
277 if (randphi) phi = Util::get_frand(0.0f,359.99999f);
278 float phitoo = params.set_default("phitoo",0.0f);
279 if ( phitoo < 0 ) throw InvalidValueException(phitoo, "Error, if you specify phitoo is must be positive");
280 Dict d;
281 d["type"] = "eman";
282 d["az"] = az;
283 d["alt"] = alt;
284 d["phi"] = phi;
285 Transform t(d);
286 v.push_back(t);
287 if ( phitoo != 0 ) {
288 if (phitoo < 0) return false;
289 else {
290 for ( float p = phitoo; p <= 360.0f-phitoo; p+= phitoo )
291 {
292 d["phi"] = fmod(phi+p,360);
293 Transform t(d);
294 v.push_back(t);
295 }
296 }
297 }
298 return true;
299}
300
301float EmanOrientationGenerator::get_az_delta(const float& delta,const float& altitude, const int) const
302{
303 // convert altitude into radians
304 float tmp = (float)(EMConsts::deg2rad * altitude);
305
306 // This is taken from EMAN1 project3d.C
307 // This wasn't working like it was supposed to. Rather than
308 // figuring it out, I'm just replacing it --steve
309/* float h=floor(360.0f/(delta*1.1547f)); // the 1.1547 makes the overall distribution more like a hexagonal mesh
310 h=(int)floor(h*sin(tmp)+.5f);
311 if (h==0) h=1;
312 h=abs(maxcsym)*floor(h/(float)abs(maxcsym)+.5f);
313 if ( h == 0 ) h = (float)maxcsym;
314 h=2.0f*M_PI/h;
315
316 return (float)(EMConsts::rad2deg*h);*/
317
318 return altitude==0?360.0f:delta/sin(tmp);
319
320}
321
322
323int EmanOrientationGenerator::get_orientations_tally(const Symmetry3D* const sym, const float& delta) const
324{
325 //FIXME THIS IS SO SIMILAR TO THE gen_orientations function that they should be probably use
326 // a common routine - SAME ISSUE FOR OTHER ORIENTATION GENERATORS
327 bool inc_mirror = params.set_default("inc_mirror",false);
328 bool breaksym = params.set_default("breaksym",false);
329 Dict delimiters = sym->get_delimiters(inc_mirror);
330 float altmax = delimiters["alt_max"];
331 float azmax = delimiters["az_max"];
332
333 float paltmin = params.set_default("alt_min",0.0f);
334 float paltmax = params.set_default("alt_max",180.0f);
335 if (altmax>paltmax) altmax=paltmax;
336
337 float alt_iterator = 0.0f;
338
339 // #If it's a h symmetry then the alt iterator starts at very close
340 // #to the altmax... the object is a h symmetry then it knows its alt_min...
341 if (sym->is_h_sym()) alt_iterator = delimiters["alt_min"];
342
343 int tally = 0;
344 while ( alt_iterator <= altmax ) {
345 float h = get_az_delta(delta,alt_iterator, sym->get_max_csym() );
346
347 // not sure what this does code taken from EMAN1 - FIXME original author add comments
348 if ( (alt_iterator > 0) && ( (azmax/h) < 2.8f) ) h = azmax / 2.1f;
349 else if (alt_iterator == 0) h = azmax;
350
351 float az_iterator = 0.0f;
352
353 float azmax_adjusted = azmax;
354 bool d_odd_mirror_flag = false;
355 get_az_max(sym,altmax, inc_mirror,alt_iterator, h,d_odd_mirror_flag, azmax_adjusted);
356 if (alt_iterator<paltmin) { alt_iterator += delta; continue; }
357
358 while ( az_iterator <= azmax_adjusted ) {
359 // FIXME: add an intelligent comment - this was copied from old code
360// if ( az_iterator > 180.0 && alt_iterator > 180.0/(2.0-0.001) && alt_iterator < 180.0/(2.0+0.001) ) {
361// az_iterator += h;
362// continue;
363// }
364
365 if (sym->is_platonic_sym()) {
366 if ( sym->is_in_asym_unit(alt_iterator, az_iterator,inc_mirror) == false ) {
367 az_iterator += h;
368 continue;
369 }
370 }
371
372 tally++;
373 if ( sym->is_h_sym() && inc_mirror && alt_iterator != (float) delimiters["alt_min"] ) {
374 tally++;
375 }
376 az_iterator += h;
377 if ( (az_iterator > azmax_adjusted) && d_odd_mirror_flag) {
378 azmax_adjusted = azmax;
379 az_iterator += azmax/2.0f;
380 }
381 }
382 alt_iterator += delta;
383 }
384
385 if (breaksym) return tally*sym->get_nsym();
386 return tally;
387}
388
389vector<Transform> EmanOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
390{
391 float delta = params.set_default("delta", 0.0f);
392 int n = params.set_default("n", 0);
393 bool breaksym = params.set_default("breaksym",false);
394 bool breaksymreal = params.set_default("breaksym_real",false);
395
396 if ( delta <= 0 && n <= 0 ) throw InvalidParameterException("Error, you must specify a positive non-zero delta or n");
397 if ( delta > 0 && n > 0 ) throw InvalidParameterException("Error, the delta and the n arguments are mutually exclusive");
398
399 if ( n > 0 ) {
400 delta = get_optimal_delta(sym,n);
401 }
402
403 bool inc_mirror = params.set_default("inc_mirror",false);
404 bool inc_mirror_real = inc_mirror;
405 if (breaksym) inc_mirror=true; // we need to enable mirror generation, then strip them out at the end, or things don't work right...
406 if (breaksymreal){
407 inc_mirror=false;
408 breaksym=true;
409 }
410 Dict delimiters = sym->get_delimiters(inc_mirror);
411 float altmax = delimiters["alt_max"];
412 float azmax = delimiters["az_max"];
413
414 float paltmin = params.set_default("alt_min",0.0f);
415 float paltmax = params.set_default("alt_max",180.0f);
416 if (altmax>paltmax) altmax=paltmax;
417
418 bool perturb = params.set_default("perturb",false); // changed to false default on 7/17/15 by steve
419
420 float alt_iterator = 0.0f;
421
422 // #If it's a h symmetry then the alt iterator starts at very close
423 // #to the altmax... the object is a h symmetry then it knows its alt_min...
424 if (sym->is_h_sym()) alt_iterator = delimiters["alt_min"];
425
426 vector<Transform> ret;
427 while ( alt_iterator <= altmax ) {
428 float h = get_az_delta(delta,alt_iterator, sym->get_max_csym() );
429
430 // not sure what this does code taken from EMAN1 - FIXME original author add comments
431 if ( (alt_iterator > 0) && ( (azmax/h) < 2.8f) ) h = azmax / 2.1f;
432 else if (alt_iterator == 0) h = azmax;
433
434 float az_iterator = 0.0f;
435
436 float azmax_adjusted = azmax;
437
438 bool d_odd_mirror_flag = false;
439 get_az_max(sym,altmax, inc_mirror,alt_iterator, h,d_odd_mirror_flag, azmax_adjusted);
440 if (alt_iterator<paltmin) { alt_iterator += delta; continue; }
441
442
443 while ( az_iterator <= azmax_adjusted ) {
444 // FIXME: add an intelligent comment - this was copied from old code
445// if ( az_iterator > 180.0 && alt_iterator > 180.0/(2.0-0.001) && alt_iterator < 180.0/(2.0+0.001) ) {
446// az_iterator += h;
447// continue;
448// }
449// // Now that I am handling the boundaries very specifically, I don't think we need
450 // the above if statement. But I am leaving it there in case I need to reconsider.
451
452 if (alt_iterator == 0 && az_iterator > 0){
453 az_iterator += h;
454 continue; // We only ever need to generate on orientation at alt=0
455 }
456
457
458 float alt_soln = alt_iterator;
459 float az_soln = az_iterator;
460
461 if (sym->is_platonic_sym()) {
462 if ( sym->is_in_asym_unit(alt_soln, az_soln,inc_mirror) == false ) {
463 az_iterator += h;
464 continue;
465 }
466 // Some objects have alignment offsets (icos and tet particularly)
467 az_soln += sym->get_az_alignment_offset();
468 }
469//printf("%f %f/n",alt_soln,az_soln);
470 if ( perturb && alt_soln != 0 ) {
471 alt_soln += Util::get_gauss_rand(0.0f,.125f*delta);
472 az_soln += Util::get_gauss_rand(0.0f,h*.125f);
473 }
474
475 add_orientation(ret,az_soln,alt_soln);
476
477 // Add helical symmetry orientations on the other side of the equator (if we're including
478 // mirror orientations)
479 if ( sym->is_h_sym() && inc_mirror && alt_iterator != (float) delimiters["alt_min"] ) {
480 add_orientation(ret, az_soln,2.0f*(float)delimiters["alt_min"]-alt_soln);
481 }
482
483 az_iterator += h;
484 if ( (az_iterator > azmax_adjusted) && d_odd_mirror_flag) {
485 azmax_adjusted = azmax;
486 az_iterator += azmax/2.0f;
487 }
488 }
489 alt_iterator += delta;
490 }
491
492 // With breaksym, values are generated for one asymmetric unit as if symmetry were imposed, then
493 // the symmetry equivalent points are generated. Used with asymmetric structures with pseudosymmetry
494 if (breaksym) {
495 // no iterators here since we are making the list longer as we go
496 int nwithsym=ret.size(); // transforms in one asym unit
497 int nsym=sym->get_nsym(); // number of asymmetric units to generate
498 for (int j=1; j<nsym; j++) {
499 Transform t=sym->get_sym(j);
500 for (int i=0; i<nwithsym; i++) {
501 ret.push_back(ret[i]*t); // add the symmetry modified transform to the end of the vector
502 }
503 }
504 if (breaksymreal) return ret;
505 // Now we get rid of anything in the bottom half of the unit sphere if requested
506 if (!inc_mirror_real) {
507 vector<Transform> ret2;
508 for (vector<Transform>::iterator t=ret.begin(); t!=ret.end(); ++t) {
509 if ((*t)[2][2]>=0) ret2.push_back(*t);
510// printf("%f\n",t[2][2]);
511 }
512 return ret2;
513 }
514 }
515
516 return ret;
517}
518
519int SingleOrientationGenerator::get_orientations_tally(const Symmetry3D* const sym, const float& delta) const
520{
521return 1;
522}
523
524vector<Transform> SingleOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
525{
526 float az = params.set_default("az", 0.0f);
527 float alt = params.set_default("alt", 0.0f);
528 float phi = params.set_default("phi", 0.0f);
529
530 float alt_iterator = 0.0f;
531
532 Transform base(Dict("type","eman","az",az,"alt",alt,"phi",phi));
533
534 vector<Transform> ret;
535 ret.push_back(base);
536
537// for (int i=0; i<sym->get_nsym(); i++) {
538// ret.push_back(base*sym->get_sym(i));
539// }
540
541
542 return ret;
543}
544
545
546vector<Transform> RandomOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
547{
548 int n = params.set_default("n", 0);
549
550 if ( n <= 0 ) throw InvalidParameterException("You must specify a positive, non zero n for the Random Orientation Generator");
551
552 bool phitoo = params.set_default("phitoo", false);
553 bool inc_mirror = params.set_default("inc_mirror", false);
554
555 vector<Transform> ret;
556
557 int i = 0;
558 Dict d("type","eman");
559 while ( i < n ){
560 float u1 = Util::get_frand(-1.0f,1.0f);
561 float u2 = Util::get_frand(-1.0f,1.0f);
562 float s = u1*u1 + u2*u2;
563 if ( s > 1.0f ) continue;
564 float alpha = 2.0f*sqrtf(1.0f-s);
565 float x = alpha * u1;
566 float y = alpha * u2;
567 float z = 2.0f*s-1.0f;
568
569 float altitude = (float)EMConsts::rad2deg*acos(z);
570 float azimuth = (float)EMConsts::rad2deg*atan2(y,x);
571
572 float phi = 0.0f;
573 if ( phitoo ) phi = Util::get_frand(0.0f,359.9999f);
574
575 d["az"] = azimuth; d["phi"] = phi; d["alt"] = altitude;
576 Transform t(d);
577
578 if ( !(sym->is_c_sym() && sym->get_nsym() == 1)) t = sym->reduce(t); //reduce doesn't make sense for C1 symmetry
579
580 if ( !sym->is_in_asym_unit(altitude,azimuth,inc_mirror) ){
581 // is_in_asym_unit has returned the wrong value!
582 // FIXME
583// cout << "warning, there is an unresolved issue - email D Woolford" << endl;
584 }
585 ret.push_back(t);
586 i++;
587 }
588 return ret;
589}
590
591int EvenOrientationGenerator::get_orientations_tally(const Symmetry3D* const sym, const float& delta) const
592{
593 bool inc_mirror = params.set_default("inc_mirror",false);
594 Dict delimiters = sym->get_delimiters(inc_mirror);
595 float altmax = delimiters["alt_max"];
596 float azmax = delimiters["az_max"];
597
598 float altmin = 0.0f;
599 // #If it's a h symmetry then the alt iterator starts at very close
600 // #to the altmax... the object is a h symmetry then it knows its alt_min...
601 if (sym->is_h_sym()) altmin = delimiters["alt_min"];
602
603 int tally = 0;
604
605 for (float alt = altmin; alt <= altmax; alt += delta) {
606 float detaz;
607 int lt;
608 if ((0.0f == alt)||(180.0f == alt)) {
609 detaz = 360.0f;
610 lt = 1;
611 } else {
612 detaz = delta/(float)sin(alt*EMConsts::deg2rad);
613 lt = int(azmax/detaz)-1;
614 if (lt < 1) lt = 1;
615 detaz = azmax/(float)lt;
616 }
617// bool d_odd_mirror_flag = false;
618// get_az_max(sym,altmax, inc_mirror,alt, lt,d_odd_mirror_flag, detaz);
619 for (int i = 0; i < lt; i++) {
620 float az = (float)i*detaz;
621 if (sym->is_platonic_sym()) {
622 if ( sym->is_in_asym_unit(alt, az,inc_mirror) == false ) continue;
623 }
624 tally++;
625 if ( sym->is_h_sym() && inc_mirror && alt != altmin ) {
626 tally++;
627 }
628 }
629 }
630
631 return tally;
632}
633
634vector<Transform> EvenOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
635{
636 float delta = params.set_default("delta", 0.0f);
637 int n = params.set_default("n", 0);
638
639 if ( delta <= 0 && n <= 0 ) throw InvalidParameterException("Error, you must specify a positive non-zero delta or n");
640 if ( delta > 0 && n > 0 ) throw InvalidParameterException("Error, the delta and the n arguments are mutually exzclusive");
641
642 if ( n > 0 ) {
643 delta = get_optimal_delta(sym,n);
644 }
645
646 bool inc_mirror = params.set_default("inc_mirror",false);
647 Dict delimiters = sym->get_delimiters(inc_mirror);
648 float altmax = delimiters["alt_max"];
649 float azmax = delimiters["az_max"];
650
651 float altmin = 0.0f;
652 // If it's a h symmetry then the alt iterator starts at very close
653 // to the altmax... the object is a h symmetry then it knows its alt_min...
654 if (sym->is_h_sym()) altmin = delimiters["alt_min"];
655
656 vector<Transform> ret;
657
658 for (float alt = altmin; alt <= altmax; alt += delta) {
659 float detaz;
660 int lt;
661 if ((0.0f == alt)||(180.0f == alt)) detaz = 360.0f;
662 else detaz = delta/(float)sin(alt*EMConsts::deg2rad);
663// bool d_odd_mirror_flag = false;
664// get_az_max(sym,altmax, inc_mirror,alt, lt,d_odd_mirror_flag, detaz);
665
666 float az = 0.0;
667 while (az<azmax) {
668 if (sym->is_platonic_sym()) {
669 if ( sym->is_in_asym_unit(alt, az, inc_mirror) ) add_orientation(ret,az+sym->get_az_alignment_offset(),alt);
670 } else add_orientation(ret,az,alt);
671 if ( sym->is_h_sym() && inc_mirror && alt != altmin ) {
672 add_orientation(ret, az, 2.0f*altmin-alt);
673 }
674 az += detaz;
675 }
676 }
677
678 return ret;
679}
680
681int SaffOrientationGenerator::get_orientations_tally(const Symmetry3D* const sym, const float& delta) const
682{
683 bool inc_mirror = params.set_default("inc_mirror",false);
684 Dict delimiters = sym->get_delimiters(inc_mirror);
685 float altmax = delimiters["alt_max"];
686 float azmax = delimiters["az_max"];
687
688 float altmin = 0.0f;
689 // #If it's a h symmetry then the alt iterator starts at very close
690 // #to the altmax... the object is a h symmetry then it knows its alt_min...
691 if (sym->is_h_sym()){
692 altmin = delimiters["alt_min"];
693 if (inc_mirror) {
694 altmin -= (float) sym->get_params()["maxtilt"];
695 }
696 }
697
698 float Deltaz = (float)(cos(altmax*EMConsts::deg2rad)-cos(altmin*EMConsts::deg2rad));
699 float s = delta*M_PI/180.0f;
700 float NFactor = 3.6f/s;
701 float wedgeFactor = fabs( Deltaz*(azmax)/720.0f) ;
702 int NumPoints = static_cast<int> (NFactor*NFactor*wedgeFactor);
703
704 int tally = 0;
705 if (!sym->is_h_sym()) ++tally;
706 float az = 0.0f;
707 float dz = (float)cos(altmin*EMConsts::deg2rad);
708 for(int i = 1; i < NumPoints; ++i ){
709 float z = dz + Deltaz* (float)i/ float(NumPoints-1);
710 float r= sqrt(1.0f-z*z);
711 az = fmod(az + delta/r,azmax);
712 float alt = (float)(acos(z)*EMConsts::rad2deg);
713 if (sym->is_platonic_sym()) {
714 if ( sym->is_in_asym_unit(alt,az,inc_mirror) == false ) continue;
715 }
716 tally++;
717 }
718
719 return tally;
720}
721
722vector<Transform> SaffOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
723{
724 float delta = params.set_default("delta", 0.0f);
725 int n = params.set_default("n", 0);
726
727 if ( delta <= 0 && n <= 0 ) throw InvalidParameterException("Error, you must specify a positive non-zero delta or n");
728 if ( delta > 0 && n > 0 ) throw InvalidParameterException("Error, the delta and the n arguments are mutually exclusive");
729
730 if ( n > 0 ) {
731 delta = get_optimal_delta(sym,n);
732 }
733
734// if ( sym->is_platonic_sym() ) return gen_platonic_orientations(sym, delta);
735
736 bool inc_mirror = params.set_default("inc_mirror",false);
737 Dict delimiters = sym->get_delimiters(inc_mirror);
738 float altmax = delimiters["alt_max"];
739 float azmax = delimiters["az_max"];
740
741 float altmin = 0.0f;
742 // #If it's a h symmetry then the alt iterator starts at very close
743 // #to the altmax... the object is a h symmetry then it knows its alt_min...
744 if (sym->is_h_sym()){
745 altmin = delimiters["alt_min"];
746 if (inc_mirror) {
747 altmin -= (float) sym->get_params()["maxtilt"];
748 }
749 }
750
751 float Deltaz = (float)(cos(altmax*EMConsts::deg2rad)-cos(altmin*EMConsts::deg2rad));
752 float s = delta*M_PI/180.0f;
753 float NFactor = 3.6f/s;
754 float wedgeFactor = fabs( Deltaz*(azmax)/720.0f) ;
755 int NumPoints = static_cast<int> (NFactor*NFactor*wedgeFactor);
756
757 vector<Transform> ret;
758
759 if (!sym->is_h_sym()) add_orientation(ret,0,0);
760 float az = 0.0f;
761 float dz = (float)cos(altmin*EMConsts::deg2rad);
762 for(int i = 1; i < NumPoints; ++i ){
763 float z = dz + Deltaz* (float)i/ float(NumPoints-1);
764 float r= sqrt(1.0f-z*z);
765 az = fmod(az + delta/r,azmax);
766 float alt = (float)(acos(z)*EMConsts::rad2deg);
767 if (sym->is_platonic_sym()) {
768 if ( sym->is_in_asym_unit(alt,az,inc_mirror) == false ) continue;
769 else {
770 az += sym->get_az_alignment_offset(); // Align to the symmetry axes
771 }
772 }
773 add_orientation(ret,az,alt);
774 }
775
776 return ret;
777}
778
779int OptimumOrientationGenerator::get_orientations_tally(const Symmetry3D* const sym, const float& delta) const
780{
781 string deltaoptname = params.set_default("use","saff");
782 Dict a;
783 a["inc_mirror"] = (bool)params.set_default("inc_mirror",false);
784 OrientationGenerator *g = Factory < OrientationGenerator >::get(deltaoptname,a);
785 if (g) {
786 int tally = g->get_orientations_tally(sym,delta);
787 delete g;
788 g = 0;
789 return tally;
790 }
791 else throw;
792}
793
794vector<Transform> OptimumOrientationGenerator::gen_orientations(const Symmetry3D* const sym) const
795{
796 float delta = params.set_default("delta", 0.0f);
797 int n = params.set_default("n", 0);
798
799 bool inc_mirror = params.set_default("inc_mirror",false);
800
801 if ( delta <= 0 && n <= 0 ) throw InvalidParameterException("Error, you must specify a positive non-zero delta or n");
802 if ( delta > 0 && n > 0 ) throw InvalidParameterException("Error, the delta and the n arguments are mutually exclusive");
803
804 string generatorname = params.set_default("use","saff");
805
806 if ( n > 0 && generatorname != RandomOrientationGenerator::NAME ) {
807 params["delta"] = get_optimal_delta(sym,n);
808 params["n"] = (int)0;
809 }
810
811 // Force the orientation generator to include the mirror - this is because
812 // We will enventually use it to generate orientations over the intire sphere
813 // which is C1 symmetry, with the inc_mirror flag set to true
814 params["inc_mirror"] = true;
815 OrientationGenerator* g = Factory < OrientationGenerator >::get(generatorname);
816 g->set_params(copy_relevant_params(g));
817
818
819 // get the starting orientation distribution
820 CSym* unit_sphere = new CSym();
821 Dict nsym; nsym["nsym"] = 1; unit_sphere->set_params(nsym);
822
823 vector<Transform> unitsphereorientations = g->gen_orientations(unit_sphere);
824 delete g; g = 0;
825 delete unit_sphere; unit_sphere = 0;
826
827 vector<Vec3f> angles = optimize_distances(unitsphereorientations);
828
829 vector<Transform> ret;
830 for (vector<Vec3f>::const_iterator it = angles.begin(); it != angles.end(); ++it ) {
831 if ( sym->is_in_asym_unit((*it)[1],(*it)[0],inc_mirror) ) {
832 add_orientation(ret,(*it)[0],(*it)[1]);
833 }
834 }
835
836 // reset the params to what they were before they were acted upon by this class
837 params["inc_mirror"] = inc_mirror;
838 params["delta"] = delta;
839 params["n"] = n;
840
841 return ret;
842}
843
844vector<Vec3f> OptimumOrientationGenerator::optimize_distances(const vector<Transform>& v) const
845{
846 vector<Vec3f> points;
847
848 for (vector<Transform>::const_iterator it = v.begin(); it != v.end(); ++it ) {
849 points.push_back(Vec3f(0,0,1)*(*it));
850 }
851
852 if ( points.size() >= 2 ) {
853 int max_it = 100;
854 float percentage = 0.01f;
855
856 for ( int i = 0; i < max_it; ++i ){
857 unsigned int p1 = 0;
858 unsigned int p2 = 1;
859
860 float distsquared = (points[p1]-points[p2]).squared_length();
861
862 // Find the nearest points
863 for(unsigned int j = 0; j < points.size(); ++j) {
864 for(unsigned int k = j+1; k < points.size(); ++k) {
865 float d = (points[j]-points[k]).squared_length();
866 if ( d < distsquared ) {
867 distsquared = d;
868 p1 = j;
869 p2 = k;
870 }
871 }
872 }
873
874 // Move them apart by a small fraction
875 Vec3f delta = percentage*(points[p2]-points[p1]);
876
877 points[p2] += delta;
878 points[p2].normalize();
879 points[p1] -= delta;
880 points[p1].normalize();
881 }
882 }
883
884 vector<Vec3f> ret;
885 for (vector<Vec3f>::const_iterator it = points.begin(); it != points.end(); ++it ) {
886 float altitude = (float)(EMConsts::rad2deg*acos((*it)[2]));
887 float azimuth = (float)(EMConsts::rad2deg*atan2((*it)[1],(*it)[0]));
888 ret.push_back(Vec3f(90.0f+azimuth,altitude,0));
889 }
890
891 return ret;
892}
893// THIS IS DWOOLFORDS FIRST SHOT AT EXTRACTING PHIL'S PLATONIC STUFF FROM SPARX
894// It didn't work.
895// vector<Transform3D> SaffOrientationGenerator::gen_platonic_orientations(const Symmetry3D* const sym, const float& delta) const
896// {
897// float scrunch = 0.9; //closeness factor to eliminate oversampling corners
898//
899// float fudge; //# fudge is a factor used to adjust phi steps
900// float m = static_cast<float>(sym->get_max_csym());
901// if ( sym->get_name() == TetrahedralSym::NAME ) fudge=0.9;
902// else if ( sym->get_name() == OctahedralSym::NAME ) fudge=0.8;
903// else if ( sym->get_name() == IcosahedralSym::NAME) fudge=0.95;
904// else throw; // this should not happen
905//
906// float n=3.0;
907// float OmegaR = 2.0*M_PI/m;
908// float cosOmega= cos(OmegaR);
909// int Edges = static_cast<int>(2.0*m*n/(2.0*(m+n)-m*n));
910// int Faces = static_cast<int>(2*Edges/n);
911// float Area = 4*M_PI/Faces/3.0; // also equals 2*pi/3 + Omega
912// float costhetac = cosOmega/(1-cosOmega);
913// float deltaRad= delta*M_PI/180.0;
914//
915// int NumPoints = static_cast<int>(Area/(deltaRad*deltaRad));
916// float fheight = 1.0f/sqrt(3.0f)/(tan(OmegaR/2.0f));
917// float z0 = costhetac; // initialize loop
918// float z = z0;
919// float phi = 0;
920// float Deltaz = (1-costhetac);
921//
922// vector<Transform3D> ret;
923// ret.push_back(Transform3D(phi, acos(z)*EMConsts::rad2deg, 0 ));
924//
925// vector<Vec3f> points;
926// points.push_back(Vec3f(sin(acos(z))*cos(phi*EMConsts::deg2rad) , sin(acos(z))*sin(phi*EMConsts::deg2rad) , z) );
927// //nLast= [ sin(acos(z))*cos(phi*piOver) , sin(acos(z))*sin(phi*piOver) , z]
928// //nVec.append(nLast)
929//
930// for(int k = 0; k < NumPoints-1; ++k ) {
931// z = z0 + Deltaz*(float)k/(float)(NumPoints-1);
932// float r= sqrt(1-z*z);
933// float phiMax =180.0*OmegaR/M_PI/2.0;
934// // Is it higher than fhat or lower
935// if (z<= fheight && false) {
936// float thetaR = acos(z);
937// float cosStuff = (cos(thetaR)/sin(thetaR))*sqrt(1. - 2 *cosOmega);
938// phiMax = 180.0*( OmegaR - acos(cosStuff))/M_PI;
939// }
940// float angleJump = fudge* delta/r;
941// phi = fmod(phi + angleJump,phiMax);
942// // anglesNew = [phi,180.0*acos(z)/pi,0.];
943// // Vec3f nNew( sin(acos(z))*cos(phi*EMConsts::deg2rad) , sin(acos(z))*sin(phi*EMConsts::deg2rad) , z);
944// // float mindiff = acos(nNew.dot(points[0]));
945// // for(unsigned int l = 0; l < points.size(); ++l ) {
946// // float dx = acos(nNew.dot(points[l]));
947// // if (dx < mindiff ) mindiff = dx;
948// // }
949// // if (mindiff > (angleJump*EMConsts::deg2rad *scrunch) ){
950// // points.push_back(nNew);
951// cout << "phi " << phi << " alt " << acos(z)*EMConsts::rad2deg << " " << z << endl;
952// ret.push_back(Transform3D(phi, acos(z)*EMConsts::rad2deg, 0 ));
953// // }
954// }
955// ret.push_back(Transform3D(0,0,0 ));
956//
957// return ret;
958// }
959
960
961
962Symmetry3D::Symmetry3D() : cached_au_planes(0),cache_size(0),num_triangles(0),au_sym_triangles() {}
963Symmetry3D::~Symmetry3D() {
964 if (cached_au_planes != 0 ) {
965 delete_au_planes();
966 }
967}
968
969void verify(const Vec3f& tmp, float * plane, const string& message )
970{
971 cout << message << " residual " << plane[0]*tmp[0]+plane[1]*tmp[1]+plane[2]*tmp[2] + plane[3] << endl;
972}
973
974Transform Symmetry3D::reduce(const Transform& t, int n) const
975{
976 // Determine which asym unit the given asym unit is in
977 int soln = in_which_asym_unit(t);
978
979 // This should never happen
980 if ( soln == -1 ) {
981 cout << "error, no solution found!" << endl;
982// throw;
983 return t;
984 }
985
986 // Get the symmetry operation corresponding to the intersection asymmetric unit
987 Transform nt = get_sym(soln);
988 // Transpose it (invert it)
989 nt.invert();
990 // Now we can transform the argument orientation into the default asymmetric unit
991 nt = t*nt;
992 // Now that we're at the default asymmetric unit, we can map into the requested asymmunit by doing this
993 if ( n != 0 ) {
994 nt = nt*get_sym(n);
995 }
996 // Done!
997 return nt;
998
999}
1000
1001int Symmetry3D::in_which_asym_unit(const Transform& t3d) const
1002{
1003 // Here it is assumed that final destination of the orientation (as encapsulated in the t3d object) is
1004 // in the z direction, so in essence we will start in the direction z and 'undo' the orientation to get the real
1005 // direction
1006 Vec3f p(0,0,1);
1007
1008 Transform o(t3d);
1009 // Orientations are alway transposed when dealing with asymmetric units, projections,etc
1010 // We take the transpose to 'undo' the transform and get the true direction of the point.
1011 o.invert();
1012 // Figure out where the point would end up. No we could just as easily not transpose and do
1013 // left multiplation (as in what occurs in the FourierReconstructor during slice insertion)
1014 p = o*p;
1015
1016 return point_in_which_asym_unit(p);
1017}
1018
1019
1020void Symmetry3D::cache_au_planes() const {
1021 if (cached_au_planes == 0 ) {
1022 vector< vector<Vec3f> > au_triangles = get_asym_unit_triangles(true);
1023 num_triangles = au_triangles.size();
1024 cache_size = get_nsym()*au_triangles.size();
1025
1026 cached_au_planes = new float*[cache_size];
1027 float** fit = cached_au_planes;
1028 for(int i =0; i < cache_size; ++i,++fit) {
1029 float *t = new float[4];
1030 *fit = t;
1031 }
1032
1033
1034 int k = 0;
1035 for(int i = 0; i < get_nsym(); ++i) {
1036
1037 for( ncit it = au_triangles.begin(); it != au_triangles.end(); ++it, ++k)
1038 {
1039 // For each given triangle
1040 vector<Vec3f> points = *it;
1041 if ( i != 0 ) {
1042 for (vector<Vec3f>::iterator iit = points.begin(); iit != points.end(); ++iit ) {
1043 // Rotate the points in the triangle so that the triangle occupies the
1044 // space of the current asymmetric unit
1045 *iit = (*iit)*get_sym(i);
1046 }
1047 }
1048
1049 au_sym_triangles.push_back(points);
1050
1051 // Determine the equation of the plane for the points, store it in plane
1052 Util::equation_of_plane(points[0],points[2],points[1],cached_au_planes[k]);
1053 }
1054 }
1055 }
1056 else throw UnexpectedBehaviorException("Attempt to generate a cache when cache exists");
1057}
1058
1059void Symmetry3D::delete_au_planes() {
1060 if (cached_au_planes == 0 ) throw UnexpectedBehaviorException("Attempt to delete a cache that does not exist");
1061 float** fit = cached_au_planes;
1062 for(int i =0; i < cache_size; ++i,++fit) {
1063 if (*fit == 0) throw UnexpectedBehaviorException("Attempt to delete a cache that does not exist");
1064 delete [] *fit;
1065 *fit = 0;
1066 }
1067
1068 delete [] cached_au_planes;
1069 cached_au_planes = 0;
1070}
1071
1072int Symmetry3D::point_in_which_asym_unit(const Vec3f& p) const
1073{
1074 if (cached_au_planes == 0) {
1075 cache_au_planes();
1076 }
1077
1078 float epsNow=0.01f;
1079 int k = 0;
1080 for(int i = 0; i < get_nsym(); ++i) {
1081 for( int j = 0; j < num_triangles; ++j,++k) {
1082 vector<Vec3f> points = au_sym_triangles[k];
1083
1084 float* plane = cached_au_planes[k];
1085 Vec3f tmp = p;
1086
1087 // Determine the intersection of p with the plane - do this by finding out how much p should be scaled by
1088 float scale = plane[0]*tmp[0]+plane[1]*tmp[1]+plane[2]*tmp[2];
1089 if ( scale != 0 )
1090 scale = -plane[3]/scale;
1091 else {
1092 // parralel!
1093 continue;
1094 }
1095
1096 // If the scale factor is less than zero, then p is definitely not in this asymmetric unit
1097 if (scale <= 0) continue;
1098
1099 // This is the intersection point
1100 Vec3f pp = tmp*scale;
1101
1102 // Now we have to see if the point p is inside the region bounded by the points, or if it is outside
1103 // If it is inside the region then p is in this asymmetric unit.
1104
1105 // This formula take from FIXME fill in once I get to work
1106 Vec3f v = points[2]-points[0];
1107 Vec3f u = points[1]-points[0];
1108 Vec3f w = pp - points[0];
1109
1110 float udotu = u.dot(u);
1111 float udotv = u.dot(v);
1112 float udotw = u.dot(w);
1113 float vdotv = v.dot(v);
1114 float vdotw = v.dot(w);
1115
1116 float d = 1.0f/(udotv*udotv - udotu*vdotv);
1117 float s = udotv*vdotw - vdotv*udotw;
1118 s *= d;
1119
1120 float t = udotv*udotw - udotu*vdotw;
1121 t *= d;
1122
1123 // We've done a few multiplications, so detect when there are tiny residuals that may throw off the final
1124 // decision
1125// printf("%d, s: %f, t: %f, d: %f\n",i,s,t,d);
1126 if (fabs(s) < Transform::ERR_LIMIT ) s = 0;
1127 if (fabs(t) < Transform::ERR_LIMIT ) t = 0;
1128
1129 // Edited by Muyuan Chen, 2015-08-19
1130 if ( fabs((fabs(s)-1.0)) < Transform::ERR_LIMIT ){
1131 if (s>0)
1132 s = 1;
1133 else
1134 s = -1;
1135 }
1136 if ( fabs((fabs(t)-1.0)) < Transform::ERR_LIMIT ){
1137 if (t>0)
1138 t = 1;
1139 else
1140 t = -1;
1141 }
1142
1143 // The final decision, if this is true then we've hit the jackpot
1144
1145
1146 if ( s >= -epsNow && t >= -epsNow && (s+t) <= 1+epsNow ) {
1147// cout << " i " << i << " j " << j << " s " << s << " t " << t << " s+t " << s+t << endl;
1148 return i;
1149 }
1150 }
1151 }
1152
1153 return -1;
1154}
1155vector<Transform> Symmetry3D::get_touching_au_transforms(bool inc_mirror) const
1156{
1157 vector<Transform> ret;
1158 vector<int> hit_cache;
1159
1160 vector<Vec3f> points = get_asym_unit_points(inc_mirror);
1161 // Warning, this is a gross hack because it is assuming that the asym_unit_points
1162 // returned by DSym are in a particular orientation with respect to symmetric axes
1163 // if the internals of DSym change it could change what we should do here...
1164 // but for the time being it will do
1165 if (inc_mirror && is_d_sym() && (get_nsym()/2 % 2 == 0)) {
1166 Dict delim = get_delimiters(false);
1167 float angle = (float)(EMConsts::deg2rad*float(delim["az_max"]));
1168 float y = -cos(angle);
1169 float x = sin(angle);
1170 points.push_back(Vec3f(x,y,0));
1171 }
1172 else if ( is_d_sym() && (get_nsym()/2 % 2 == 1)) {
1173 Dict delim = get_delimiters(false);
1174 float angle = float(delim["az_max"])/2.0f;
1175// cout << "Odd dsym using " << angle << endl;
1176 angle *= (float)EMConsts::deg2rad;
1177 float y = -cos(angle);
1178 float x = sin(angle);
1179 points.push_back(Vec3f(x,y,0));
1180
1181 if ( inc_mirror ) {
1182 angle = 3.0f*(float(delim["az_max"]))/2.0f;
1183 angle *= (float)EMConsts::deg2rad;
1184 float y = -cos(angle);
1185 float x = sin(angle);
1186 points.push_back(Vec3f(x,y,0));
1187 }
1188 }
1189
1190 typedef vector<Vec3f>::const_iterator const_point_it;
1191 for(const_point_it point = points.begin(); point != points.end(); ++point ) {
1192
1193 for(int i = 1; i < get_nsym(); ++i) {
1194
1195 if ( find(hit_cache.begin(),hit_cache.end(),i) != hit_cache.end() ) continue;
1196 Transform t = get_sym(i);
1197 Vec3f result = (*point)*t;
1198
1199 if (is_platonic_sym()) {
1200 for(const_point_it tmp = points.begin(); tmp != points.end(); ++tmp ) {
1201 Vec3f tt = result-(*tmp);
1202 if (tt.squared_length() < 0.01f) {
1203 hit_cache.push_back(i);
1204 ret.push_back(t);
1205 }
1206
1207 }
1208 }
1209 else {
1210 result -= *point;
1211 if (result.squared_length() < 0.05f) {
1212 hit_cache.push_back(i);
1213 ret.push_back(t);
1214 }
1215 }
1216 }
1217
1218 }
1219
1220 return ret;
1221}
1222
1223
1224vector<Transform> Symmetry3D::get_syms() const
1225{
1226
1227 vector<Transform> ret;
1228// if (t.is_identity()) {
1229 for(int i = 0; i < get_nsym(); ++i ) {
1230 ret.push_back(get_sym(i));
1231 }
1232// } else {
1233// for(int i = 0; i < get_nsym(); ++i ) {
1234// ret.push_back(get_sym(i)*t);
1235// }
1236// }
1237 return ret;
1238}
1239
1240vector<Transform> Symmetry3D::get_symmetries(const string& symmetry)
1241{
1242 Symmetry3D* sym = Factory<Symmetry3D>::get(Util::str_to_lower(symmetry));
1243 vector<Transform> ret = sym->get_syms();
1244 delete sym;
1245 return ret;
1246}
1247
1248// C Symmetry stuff
1249Dict CSym::get_delimiters(const bool inc_mirror) const {
1250 Dict returnDict;
1251 // Get the parameters of interest
1252 int nsym = params.set_default("nsym",0);
1253 if ( nsym <= 0 ) throw InvalidValueException(nsym,"Error, you must specify a positive non zero nsym");
1254
1255 if ( inc_mirror ) returnDict["alt_max"] = 180.0f;
1256 else returnDict["alt_max"] = 90.0f;
1257
1258 returnDict["az_max"] = 360.0f/(float)nsym;
1259
1260 return returnDict;
1261}
1262
1263bool CSym::is_in_asym_unit(const float& altitude, const float& azimuth, const bool inc_mirror = false) const
1264{
1265 Dict d = get_delimiters(inc_mirror);
1266 float alt_max = d["alt_max"];
1267 float az_max = d["az_max"];
1268
1269 int nsym = params.set_default("nsym",0);
1270 if ( nsym != 1 && azimuth < 0) return false;
1271 if ( altitude <= alt_max && azimuth < az_max ) return true;
1272 return false;
1273}
1274
1275vector<vector<Vec3f> > CSym::get_asym_unit_triangles(bool inc_mirror ) const{
1276 vector<Vec3f> v = get_asym_unit_points(inc_mirror);
1277 int nsym = params.set_default("nsym",0);
1278
1279 vector<vector<Vec3f> > ret;
1280 if (v.size() == 0) return ret; // nsym == 1 and inc_mirror == true, this is the entire sphere!
1281 if (nsym == 1 && !inc_mirror) {
1282 Vec3f z(0,0,1);
1283 vector<Vec3f> tmp;
1284 tmp.push_back(z);
1285 tmp.push_back(v[1]);
1286 tmp.push_back(v[0]);
1287 ret.push_back(tmp);
1288
1289 vector<Vec3f> tmp2;
1290 tmp2.push_back(z);
1291 tmp2.push_back(v[2]);
1292 tmp2.push_back(v[1]);
1293 ret.push_back(tmp2);
1294
1295 vector<Vec3f> tmp3;
1296 tmp3.push_back(z);
1297 tmp3.push_back(v[3]);
1298 tmp3.push_back(v[2]);
1299 ret.push_back(tmp3);
1300
1301 vector<Vec3f> tmp4;
1302 tmp4.push_back(z);
1303 tmp4.push_back(v[0]);
1304 tmp4.push_back(v[3]);
1305 ret.push_back(tmp4);
1306 }
1307 else if (nsym == 2 && inc_mirror) {
1308 Vec3f x(1,0,0);
1309 vector<Vec3f> tmp;
1310 tmp.push_back(v[1]);
1311 tmp.push_back(v[0]);
1312 tmp.push_back(x);
1313 ret.push_back(tmp);
1314
1315 vector<Vec3f> tmp2;
1316 tmp2.push_back(v[2]);
1317 tmp2.push_back(v[1]);
1318 tmp2.push_back(x);
1319 ret.push_back(tmp2);
1320
1321 vector<Vec3f> tmp3;
1322 tmp3.push_back(v[3]);
1323 tmp3.push_back(v[2]);
1324 tmp3.push_back(x);
1325 ret.push_back(tmp3);
1326
1327 vector<Vec3f> tmp4;
1328 tmp4.push_back(v[0]);
1329 tmp4.push_back(v[3]);
1330 tmp4.push_back(x);
1331 ret.push_back(tmp4);
1332 }
1333 else if (nsym == 2 && !inc_mirror) {
1334 vector<Vec3f> tmp;
1335 tmp.push_back(v[0]);
1336 tmp.push_back(v[2]);
1337 tmp.push_back(v[1]);
1338 ret.push_back(tmp);
1339
1340 vector<Vec3f> tmp2;
1341 tmp2.push_back(v[2]);
1342 tmp2.push_back(v[0]);
1343 tmp2.push_back(v[3]);
1344 ret.push_back(tmp2);
1345 }
1346 else if (v.size() == 3) {
1347 vector<Vec3f> tmp;
1348 tmp.push_back(v[0]);
1349 tmp.push_back(v[2]);
1350 tmp.push_back(v[1]);
1351 ret.push_back(tmp);
1352 }
1353 else if (v.size() == 4) {
1354 vector<Vec3f> tmp;
1355 tmp.push_back(v[0]);
1356 tmp.push_back(v[3]);
1357 tmp.push_back(v[1]);
1358 ret.push_back(tmp);
1359
1360 vector<Vec3f> tmp2;
1361 tmp2.push_back(v[1]);
1362 tmp2.push_back(v[3]);
1363 tmp2.push_back(v[2]);
1364 ret.push_back(tmp2);
1365 }
1366
1367 return ret;
1368}
1369
1370vector<Vec3f> CSym::get_asym_unit_points(bool inc_mirror) const
1371{
1372 Dict delim = get_delimiters(inc_mirror);
1373 int nsym = params.set_default("nsym",0);
1374 vector<Vec3f> ret;
1375
1376 if ( nsym == 1 ) {
1377 if (inc_mirror == false ) {
1378 ret.push_back(Vec3f(0,-1,0));
1379 ret.push_back(Vec3f(1,0,0));
1380 ret.push_back(Vec3f(0,1,0));
1381 ret.push_back(Vec3f(-1,0,0));
1382 }
1383 // else return ret; // an empty vector! this is fine
1384 }
1385 else if (nsym == 2 && !inc_mirror) {
1386 ret.push_back(Vec3f(0,0,1));
1387 ret.push_back(Vec3f(0,-1,0));
1388 ret.push_back(Vec3f(1,0,0));
1389 ret.push_back(Vec3f(0,1,0));
1390 }
1391 else {
1392 ret.push_back(Vec3f(0,0,1));
1393 ret.push_back(Vec3f(0,-1,0));
1394 if (inc_mirror == true) {
1395 ret.push_back(Vec3f(0,0,-1));
1396 }
1397 float angle = (float)(EMConsts::deg2rad*float(delim["az_max"]));
1398 float y = -cos(angle);
1399 float x = sin(angle);
1400 ret.push_back(Vec3f(x,y,0));
1401 }
1402
1403 return ret;
1404
1405}
1406
1407Transform CSym::get_sym(const int n) const {
1408 int nsym = params.set_default("nsym",0);
1409 if ( nsym <= 0 ) throw InvalidValueException(n,"Error, you must specify a positive non zero nsym");
1410
1411 Dict d("type","eman");
1412 // courtesy of Phil Baldwin
1413 d["az"] = (n%nsym) * 360.0f / nsym;
1414 d["alt"] = 0.0f;
1415 d["phi"] = 0.0f;
1416 return Transform(d);
1417}
1418
1419// D symmetry stuff
1420Dict DSym::get_delimiters(const bool inc_mirror) const {
1421 Dict returnDict;
1422
1423 // Get the parameters of interest
1424 int nsym = params.set_default("nsym",0);
1425 if ( nsym <= 0 ) throw InvalidValueException(nsym,"Error, you must specify a positive non zero nsym");
1426
1427 returnDict["alt_max"] = 90.0f;
1428
1429 if ( inc_mirror ) returnDict["az_max"] = 360.0f/(float)nsym;
1430 else returnDict["az_max"] = 180.0f/(float)nsym;
1431
1432 return returnDict;
1433}
1434
1435bool DSym::is_in_asym_unit(const float& altitude, const float& azimuth, const bool inc_mirror = false) const
1436{
1437 Dict d = get_delimiters(inc_mirror);
1438 float alt_max = d["alt_max"];
1439 float az_max = d["az_max"];
1440
1441 int nsym = params.set_default("nsym",0);
1442
1443 if ( nsym == 1 && inc_mirror ) {
1444 if (altitude >= 0 && altitude <= alt_max && azimuth < az_max ) return true;
1445 }
1446 else {
1447 if ( altitude >= 0 && altitude <= alt_max && azimuth < az_max && azimuth >= 0 ) return true;
1448 }
1449 return false;
1450}
1451
1452Transform DSym::get_sym(const int n) const
1453{
1454 int nsym = 2*params.set_default("nsym",0);
1455 if ( nsym <= 0 ) throw InvalidValueException(n,"Error, you must specify a positive non zero nsym");
1456
1457 Dict d("type","eman");
1458 // courtesy of Phil Baldwin
1459 if (n >= nsym / 2) {
1460 d["az"] = ( (n%nsym) - nsym/2) * 360.0f / (nsym / 2);
1461 d["alt"] = 180.0f;
1462 d["phi"] = 0.0f;
1463 }
1464 else {
1465 d["az"] = (n%nsym) * 360.0f / (nsym / 2);
1466 d["alt"] = 0.0f;
1467 d["phi"] = 0.0f;
1468 }
1469 return Transform(d);
1470}
1471
1472vector<vector<Vec3f> > DSym::get_asym_unit_triangles(bool inc_mirror) const{
1473 vector<Vec3f> v = get_asym_unit_points(inc_mirror);
1474 int nsym = params.set_default("nsym",0);
1475 vector<vector<Vec3f> > ret;
1476 if ( (nsym == 1 && inc_mirror == false) || (nsym == 2 && inc_mirror)) {
1477 vector<Vec3f> tmp;
1478 tmp.push_back(v[0]);
1479 tmp.push_back(v[2]);
1480 tmp.push_back(v[1]);
1481 ret.push_back(tmp);
1482
1483 vector<Vec3f> tmp2;
1484 tmp2.push_back(v[2]);
1485 tmp2.push_back(v[0]);
1486 tmp2.push_back(v[3]);
1487 ret.push_back(tmp2);
1488 }
1489 else if (nsym == 1) {
1490 Vec3f z(0,0,1);
1491 vector<Vec3f> tmp;
1492 tmp.push_back(z);
1493 tmp.push_back(v[1]);
1494 tmp.push_back(v[0]);
1495 ret.push_back(tmp);
1496
1497 vector<Vec3f> tmp2;
1498 tmp2.push_back(z);
1499 tmp2.push_back(v[2]);
1500 tmp2.push_back(v[1]);
1501 ret.push_back(tmp2);
1502
1503 vector<Vec3f> tmp3;
1504 tmp3.push_back(z);
1505 tmp3.push_back(v[3]);
1506 tmp3.push_back(v[2]);
1507 ret.push_back(tmp3);
1508
1509 vector<Vec3f> tmp4;
1510 tmp4.push_back(z);
1511 tmp4.push_back(v[0]);
1512 tmp4.push_back(v[3]);
1513 ret.push_back(tmp4);
1514 }
1515 else {
1516// if v.size() == 3
1517 vector<Vec3f> tmp;
1518 tmp.push_back(v[0]);
1519 tmp.push_back(v[2]);
1520 tmp.push_back(v[1]);
1521 ret.push_back(tmp);
1522 }
1523
1524 return ret;
1525}
1526
1527vector<Vec3f> DSym::get_asym_unit_points(bool inc_mirror) const
1528{
1529 Dict delim = get_delimiters(inc_mirror);
1530
1531 vector<Vec3f> ret;
1532 int nsym = params.set_default("nsym",0);
1533 if ( nsym == 1 ) {
1534 if (inc_mirror == false ) {
1535 ret.push_back(Vec3f(0,0,1));
1536 ret.push_back(Vec3f(0,-1,0));
1537 ret.push_back(Vec3f(1,0,0));
1538 ret.push_back(Vec3f(0,1,0));
1539 }
1540 else {
1541 ret.push_back(Vec3f(0,-1,0));
1542 ret.push_back(Vec3f(1,0,0));
1543 ret.push_back(Vec3f(0,1,0));
1544 ret.push_back(Vec3f(-1,0,0));
1545 }
1546 }
1547 else if ( nsym == 2 && inc_mirror ) {
1548 ret.push_back(Vec3f(0,0,1));
1549 ret.push_back(Vec3f(0,-1,0));
1550 ret.push_back(Vec3f(1,0,0));
1551 ret.push_back(Vec3f(0,1,0));
1552 }
1553 else {
1554 float angle = (float)(EMConsts::deg2rad*float(delim["az_max"]));
1555 ret.push_back(Vec3f(0,0,1));
1556 ret.push_back(Vec3f(0,-1,0));
1557 float y = -cos(angle);
1558 float x = sin(angle);
1559 ret.push_back(Vec3f(x,y,0));
1560 }
1561
1562 return ret;
1563
1564}
1565
1566// H symmetry stuff
1567Dict HSym::get_delimiters(const bool) const {
1568 Dict returnDict;
1569
1570 // Get the parameters of interest
1571 int nsym = params.set_default("nsym",0);
1572 if ( nsym <= 0 ) throw InvalidValueException(nsym,"Error, you must specify a positive non zero nsym");
1573
1574 float maxtilt = params.set_default("maxtilt",90.0f);
1575
1576 returnDict["alt_max"] = 90.0f;
1577
1578 returnDict["alt_min"] = 90.0f - maxtilt;
1579
1580 returnDict["az_max"] = 360.0f;
1581
1582 return returnDict;
1583}
1584
1585bool HSym::is_in_asym_unit(const float& altitude, const float& azimuth, const bool inc_mirror = false) const
1586{
1587 Dict d = get_delimiters(inc_mirror);
1588 float alt_max = d["alt_max"];
1589 float alt_min = d["alt_min"];
1590
1591 if (inc_mirror) {
1592 float e = params.set_default("maxtilt",90.0f);
1593 alt_min -= e;
1594 }
1595
1596 float az_max = d["az_max"];
1597
1598 if ( altitude >=alt_min && altitude <= alt_max && azimuth <= az_max && azimuth >= 0 ) return true;
1599 return false;
1600}
1601
1602vector<vector<Vec3f> > HSym::get_asym_unit_triangles(bool ) const{
1603
1604 vector<vector<Vec3f> > ret;
1605 return ret;
1606}
1607
1608vector<Vec3f> HSym::get_asym_unit_points(bool inc_mirror) const
1609{
1610 vector<Vec3f> ret;
1611
1612 Dict delim = get_delimiters(inc_mirror);
1613 int nsym = params.set_default("nsym",1);
1614 float az = -(float)delim["az_max"];
1615
1616
1617 bool tracing_arcs = false;
1618
1619
1620 if ( !tracing_arcs) {
1621 Vec3f a(0,-1,0);
1622 ret.push_back(a);
1623
1624 if ( nsym > 2 ) {
1625 Dict d("type","eman");
1626 d["phi"] = 0.0f;
1627 d["alt"] = 0.0f;
1628 d["az"] = az;
1629 Vec3f b = Transform(d)*a;
1630 ret.push_back(b);
1631 }
1632 else
1633 {
1634 ret.push_back(Vec3f(1,0,0));
1635
1636 ret.push_back(Vec3f(0,1,0));
1637
1638 if ( nsym == 1 ) {
1639 ret.push_back(Vec3f(-1,0,0));
1640 ret.push_back(a);
1641 }
1642 }
1643 }
1644 return ret;
1645
1646}
1647
1648Transform HSym::get_sym(const int n) const
1649{
1650 int nstart=params["nstart"];
1651 //int nsym=params["nsym"];
1652 float apix = params.set_default("apix",1.0f);
1653 float daz= params["daz"];
1654 float tz=params["tz"];
1655 float dz=tz/apix;
1656 Dict d("type","eman");
1657
1658 // courtesy of Phil Baldwin
1659 //d["az"] = (n%nsym) * 360.0f / nsym;
1660 //d["az"]=(((int) n/hsym)%nstart)*360.f/nstart+(n%hsym)*daz;
1661 //d["az"] = n * daz;
1662 int ii=(n+1)/2;
1663 if (n>1 && n%2==0) ii*=-1; // extend to both directions
1664 d["az"]=(ii%nstart)*(360.0/nstart)+floor(float(ii)/nstart)*daz; // corrected by steve, 7/21/11. No dependency on nsym
1665 d["alt"] = 0.0f;
1666 d["phi"] = 0.0f;
1667 Transform ret(d);
1668 ret.set_trans(0,0,(ii/nstart)*dz);
1669 return ret;
1670}
1671
1672// Generic platonic symmetry stuff
1673void PlatonicSym::init()
1674{
1675 //See the manuscript "The Transform Class in Sparx and EMAN2", Baldwin & Penczek 2007. J. Struct. Biol. 157 (250-261)
1676 //In particular see pages 257-259
1677 //cap_sig is capital sigma in the Baldwin paper
1678 float cap_sig = 2.0f*M_PI/ get_max_csym();
1679 //In EMAN2 projection cap_sig is really az_max
1680 platonic_params["az_max"] = cap_sig;
1681
1682 // Alpha is the angle between (immediately) neighborhing 3 fold axes of symmetry
1683 // This follows the conventions in the Baldwin paper
1684 float alpha = acos(1.0f/(sqrtf(3.0f)*tan(cap_sig/2.0f)));
1685 // In EMAN2 projection alpha is really al_maz
1686 platonic_params["alt_max"] = alpha;
1687
1688 // This is half of "theta_c" as in the conventions of the Balwin paper. See also http://blake.bcm.edu/emanwiki/EMAN2/Symmetry.
1689 platonic_params["theta_c_on_two"] = 1.0f/2.0f*acos( cos(cap_sig)/(1.0f-cos(cap_sig)));
1690
1691}
1692
1693
1694Dict PlatonicSym::get_delimiters(const bool inc_mirror) const
1695{
1696 Dict ret;
1697 ret["az_max"] = EMConsts::rad2deg * (float) platonic_params["az_max"];
1698 // For icos and oct symmetries, excluding the mirror means halving az_maz
1699 if ( inc_mirror == false )
1700 if ( get_name() == IcosahedralSym::NAME || get_name() == OctahedralSym::NAME )
1701 ret["az_max"] = 0.5f*EMConsts::rad2deg * (float) platonic_params["az_max"];
1702 //else
1703 //the alt_max variable should probably be altered if the symmetry is tet, but
1704 //this is taken care of in TetSym::is_in_asym_unit
1705
1706 ret["alt_max"] = (float)(EMConsts::rad2deg * (float) platonic_params["alt_max"]);
1707 return ret;
1708}
1709
1710//.Warning, this function only returns valid answers for octahedral and icosahedral symmetries.
1711bool PlatonicSym::is_in_asym_unit(const float& altitude, const float& azimuth, const bool inc_mirror) const
1712{
1713 Dict d = get_delimiters(inc_mirror);
1714 float alt_max = d["alt_max"];
1715 float az_max = d["az_max"];
1716
1717 if ( altitude >= 0 && altitude <= alt_max && azimuth < az_max && azimuth >= 0) {
1718
1719 // Convert azimuth to radians
1720 float tmpaz = (float)(EMConsts::deg2rad * azimuth);
1721
1722 float cap_sig = platonic_params["az_max"];
1723 float alt_max = platonic_params["alt_max"];
1724 tmpaz = Util::get_min( tmpaz, cap_sig-tmpaz);
1725
1726 // convert altitude to radians
1727 float tmpalt = (float)(EMConsts::deg2rad * altitude);
1728 if ( platonic_alt_lower_bound(tmpaz, alt_max ) > tmpalt ) return true;
1729 /*if ( inc_mirror == false ) {
1730 if ( tmpaz > cap_sig/2.0f) return false; // this if is always false so I removed it PAP
1731 else return true;
1732 } else return true;*/
1733 else return false;
1734 } else return false;
1735}
1736
1737float PlatonicSym::platonic_alt_lower_bound(const float& azimuth, const float& alpha) const
1738{
1739 float cap_sig = platonic_params["az_max"];
1740 float theta_c_on_two = platonic_params["theta_c_on_two"];
1741
1742 float baldwin_lower_alt_bound = sin(cap_sig/2.0f-azimuth)/tan(theta_c_on_two);
1743 baldwin_lower_alt_bound += sin(azimuth)/tan(alpha);
1744 baldwin_lower_alt_bound *= 1/sin(cap_sig/2.0f);
1745 baldwin_lower_alt_bound = atan(1/baldwin_lower_alt_bound);
1746
1747 return baldwin_lower_alt_bound;
1748}
1749
1750vector<vector<Vec3f> > PlatonicSym::get_asym_unit_triangles(bool inc_mirror) const{
1751 vector<Vec3f> v = get_asym_unit_points(inc_mirror);
1752 vector<vector<Vec3f> > ret;
1753 if (v.size() == 3) {
1754 vector<Vec3f> tmp;
1755 tmp.push_back(v[0]);
1756 tmp.push_back(v[2]);
1757 tmp.push_back(v[1]);
1758 ret.push_back(tmp);
1759 }
1760 else /* v.size() == 4*/ {
1761 vector<Vec3f> tmp;
1762 tmp.push_back(v[0]);
1763 tmp.push_back(v[2]);
1764 tmp.push_back(v[1]);
1765 ret.push_back(tmp);
1766
1767 vector<Vec3f> tmp2;
1768 tmp2.push_back(v[0]);
1769 tmp2.push_back(v[3]);
1770 tmp2.push_back(v[2]);
1771 ret.push_back(tmp2);
1772 }
1773
1774 return ret;
1775}
1776
1777vector<Vec3f> PlatonicSym::get_asym_unit_points(bool inc_mirror) const
1778{
1779 vector<Vec3f> ret;
1780
1781 Vec3f b = Vec3f(0,0,1);
1782 ret.push_back(b);
1783 float theta_c_on_two = (float)platonic_params["theta_c_on_two"]; // already in radians
1784 float theta_c = 2*theta_c_on_two;
1785
1786 Vec3f c_on_two = Vec3f(0,-sin(theta_c_on_two),cos(theta_c_on_two));
1787 Vec3f c = Vec3f(0,-sin(theta_c),cos(theta_c));
1788 ret.push_back(c_on_two);
1789
1790 float cap_sig = platonic_params["az_max"];
1791 Vec3f a = Vec3f(sin(theta_c)*sin(cap_sig),-sin(theta_c)*cos(cap_sig),cos(theta_c));
1792
1793 Vec3f f = a+b+c;
1794 f.normalize();
1795
1796 ret.push_back(f);
1797
1798 if ( inc_mirror ) {
1799 Vec3f a_on_two = Vec3f(sin(theta_c_on_two)*sin(cap_sig),-sin(theta_c_on_two)*cos(cap_sig),cos(theta_c_on_two));
1800 ret.push_back(a_on_two);
1801 }
1802
1803 if ( get_az_alignment_offset() != 0 ) {
1804 Dict d("type","eman");
1805 d["az"] = get_az_alignment_offset();
1806 d["phi"] = 0.0f;
1807 d["alt"] = 0.0f;
1808 Transform t(d);
1809 for (vector<Vec3f>::iterator it = ret.begin(); it != ret.end(); ++it )
1810 {
1811 *it = (*it)*t;
1812 }
1813 }
1814 //
1815 return ret;
1816
1817}
1818
1819float IcosahedralSym::get_az_alignment_offset() const { return 234.0; } // This offset positions a 3 fold axis on the positive x axis
1820
1821Transform IcosahedralSym::get_sym(const int n) const
1822{
1823 // These rotations courtesy of Phil Baldwin
1824 float lvl0 = 0.0f; // there is one pentagon on top; five-fold along z
1825 float lvl1 = EMConsts::rad2deg*atan(2.0f);// 63.4349; // that is atan(2) // there are 5 pentagons with centers at this height (angle)
1826 float lvl2 = 180.0f - lvl1;//116.5651; //that is 180-lvl1 // there are 5 pentagons with centers at this height (angle)
1827 float lvl3 = 180.0f;
1828
1829 float ICOS[180] = { // This is with a pentagon normal to z
1830 0,lvl0,0, 0,lvl0,288, 0,lvl0,216, 0,lvl0,144, 0,lvl0,72,
1831 0,lvl1,36, 0,lvl1,324, 0,lvl1,252, 0,lvl1,180, 0,lvl1,108,
1832 72,lvl1,36, 72,lvl1,324, 72,lvl1,252, 72,lvl1,180, 72,lvl1,108,
1833 144,lvl1,36, 144,lvl1,324, 144,lvl1,252, 144,lvl1,180, 144,lvl1,108,
1834 216,lvl1,36, 216,lvl1,324, 216,lvl1,252, 216,lvl1,180, 216,lvl1,108,
1835 288,lvl1,36, 288,lvl1,324, 288,lvl1,252, 288,lvl1,180, 288,lvl1,108,
1836 36,lvl2,0, 36,lvl2,288, 36,lvl2,216, 36,lvl2,144, 36,lvl2,72,
1837 108,lvl2,0, 108,lvl2,288, 108,lvl2,216, 108,lvl2,144, 108,lvl2,72,
1838 180,lvl2,0, 180,lvl2,288, 180,lvl2,216, 180,lvl2,144, 180,lvl2,72,
1839 252,lvl2,0, 252,lvl2,288, 252,lvl2,216, 252,lvl2,144, 252,lvl2,72,
1840 324,lvl2,0, 324,lvl2,288, 324,lvl2,216, 324,lvl2,144, 324,lvl2,72,
1841 0,lvl3,0, 0,lvl3,288, 0,lvl3,216, 0,lvl3,144, 0,lvl3,72
1842 };
1843
1844 int idx = n % 60;
1845 Dict d("type","eman");
1846// Transform3D ret;
1847 if (get_az_alignment_offset() == 234.0) {
1848 d["az"] = ICOS[idx * 3 ]+90;
1849 d["alt"] = ICOS[idx * 3 + 1];
1850 d["phi"] = ICOS[idx * 3 + 2]-90;
1851// ret.set_rotation((float)ICOS[idx * 3 ]+90,(float)ICOS[idx * 3 + 1], (float)ICOS[idx * 3 + 2]-90);
1852 }
1853 else {
1854 d["az"] = ICOS[idx * 3 ];
1855 d["alt"] = ICOS[idx * 3 + 1];
1856 d["phi"] = ICOS[idx * 3 + 2];
1857// ret.set_rotation((float)ICOS[idx * 3 ],(float)ICOS[idx * 3 + 1], (float)ICOS[idx * 3 + 2]);
1858 }
1859
1860// ret.set_rotation((float)ICOS[idx * 3 ],(float)ICOS[idx * 3 + 1], (float)ICOS[idx * 3 + 2]);
1861// if ( get_az_alignment_offset() != 0 ) {
1862// Transform3D t(get_az_alignment_offset(),0,0);
1863// ret = t*ret;
1864// }
1865 return Transform(d);
1866
1867}
1868
1869float Icosahedral2Sym::get_az_alignment_offset() const { return 234.0; } // This offset positions a 3 fold axis on the positive x axis (??? copied from IcosahedralSym)
1870
1871Transform Icosahedral2Sym::get_sym(const int n) const
1872{
1873 static float matrices[60*9] = {
1874 1, 0, 0, 0, 1, 0, 0, 0, 1,
1875 0.30902, -0.80902, 0.5, 0.80902, 0.5, 0.30902, -0.5, 0.30902, 0.80902,
1876 -0.80902, -0.5, 0.30902, 0.5, -0.30902, 0.80902, -0.30902, 0.80902, 0.5,
1877 -0.80902, 0.5, -0.30902, -0.5, -0.30902, 0.80902, 0.30902, 0.80902, 0.5,
1878 0.30902, 0.80902, -0.5, -0.80902, 0.5, 0.30902, 0.5, 0.30902, 0.80902,
1879 -1, 0, 0, 0, -1, 0, 0, 0, 1,
1880 -0.30902, -0.80902, 0.5, 0.80902, -0.5, -0.30902, 0.5, 0.30902, 0.80902,
1881 0.30902, 0.80902, 0.5, -0.80902, 0.5, -0.30902, -0.5, -0.30902, 0.80902,
1882 -0.30902, 0.80902, 0.5, -0.80902, -0.5, 0.30902, 0.5, -0.30902, 0.80902,
1883 -0.5, 0.30902, 0.80902, 0.30902, -0.80902, 0.5, 0.80902, 0.5, 0.30902,
1884 0.5, -0.30902, 0.80902, -0.30902, 0.80902, 0.5, -0.80902, -0.5, 0.30902,
1885 0.80902, 0.5, 0.30902, -0.5, 0.30902, 0.80902, 0.30902, -0.80902, 0.5,
1886 0.80902, -0.5, -0.30902, 0.5, 0.30902, 0.80902, -0.30902, -0.80902, 0.5,
1887 0.5, 0.30902, -0.80902, 0.30902, 0.80902, 0.5, 0.80902, -0.5, 0.30902,
1888 -0.5, -0.30902, -0.80902, -0.30902, -0.80902, 0.5, -0.80902, 0.5, 0.30902,
1889 -0.30902, -0.80902, -0.5, 0.80902, -0.5, 0.30902, -0.5, -0.30902, 0.80902,
1890 0.30902, -0.80902, -0.5, 0.80902, 0.5, -0.30902, 0.5, -0.30902, 0.80902,
1891 -0.30902, 0.80902, -0.5, -0.80902, -0.5, -0.30902, -0.5, 0.30902, 0.80902,
1892 0.5, -0.30902, -0.80902, -0.30902, 0.80902, -0.5, 0.80902, 0.5, 0.30902,
1893 -0.5, 0.30902, -0.80902, 0.30902, -0.80902, -0.5, -0.80902, -0.5, 0.30902,
1894 -0.80902, -0.5, -0.30902, 0.5, -0.30902, -0.80902, 0.30902, -0.80902, 0.5,
1895 0.80902, 0.5, -0.30902, -0.5, 0.30902, -0.80902, -0.30902, 0.80902, 0.5,
1896 0.80902, -0.5, 0.30902, 0.5, 0.30902, -0.80902, 0.30902, 0.80902, 0.5,
1897 -0.80902, 0.5, 0.30902, -0.5, -0.30902, -0.80902, -0.30902, -0.80902, 0.5,
1898 -0.5, -0.30902, 0.80902, -0.30902, -0.80902, -0.5, 0.80902, -0.5, 0.30902,
1899 0.5, 0.30902, 0.80902, 0.30902, 0.80902, -0.5, -0.80902, 0.5, 0.30902,
1900 0, 0, 1, 1, 0, 0, 0, 1, 0,
1901 0, 1, 0, 0, 0, 1, 1, 0, 0,
1902 0, -1, 0, 0, 0, 1, -1, 0, 0,
1903 0, 0, -1, -1, 0, 0, 0, 1, 0,
1904 0, -1, 0, 0, 0, -1, 1, 0, 0,
1905 0, 1, 0, 0, 0, -1, -1, 0, 0,
1906 -0.80902, -0.5, 0.30902, -0.5, 0.30902, -0.80902, 0.30902, -0.80902, -0.5,
1907 0.80902, -0.5, -0.30902, -0.5, -0.30902, -0.80902, 0.30902, 0.80902, -0.5,
1908 0.5, 0.30902, -0.80902, -0.30902, -0.80902, -0.5, -0.80902, 0.5, -0.30902,
1909 -0.30902, -0.80902, -0.5, -0.80902, 0.5, -0.30902, 0.5, 0.30902, -0.80902,
1910 -0.80902, 0.5, -0.30902, 0.5, 0.30902, -0.80902, -0.30902, -0.80902, -0.5,
1911 -0.5, -0.30902, -0.80902, 0.30902, 0.80902, -0.5, 0.80902, -0.5, -0.30902,
1912 -0.5, 0.30902, -0.80902, -0.30902, 0.80902, 0.5, 0.80902, 0.5, -0.30902,
1913 0, 0, -1, 1, 0, 0, 0, -1, 0,
1914 -0.80902, 0.5, 0.30902, 0.5, 0.30902, 0.80902, 0.30902, 0.80902, -0.5,
1915 0.80902, 0.5, -0.30902, 0.5, -0.30902, 0.80902, 0.30902, -0.80902, -0.5,
1916 -0.30902, 0.80902, -0.5, 0.80902, 0.5, 0.30902, 0.5, -0.30902, -0.80902,
1917 0.5, -0.30902, -0.80902, 0.30902, -0.80902, 0.5, -0.80902, -0.5, -0.30902,
1918 -0.80902, -0.5, -0.30902, -0.5, 0.30902, 0.80902, -0.30902, 0.80902, -0.5,
1919 -0.30902, -0.80902, 0.5, -0.80902, 0.5, 0.30902, -0.5, -0.30902, -0.80902,
1920 -0.30902, 0.80902, 0.5, 0.80902, 0.5, -0.30902, -0.5, 0.30902, -0.80902,
1921 1, 0, 0, 0, -1, 0, 0, 0, -1,
1922 0.30902, 0.80902, -0.5, 0.80902, -0.5, -0.30902, -0.5, -0.30902, -0.80902,
1923 0.30902, -0.80902, -0.5, -0.80902, -0.5, 0.30902, -0.5, 0.30902, -0.80902,
1924 -1, 0, 0, 0, 1, 0, 0, 0, -1,
1925 0.80902, 0.5, 0.30902, 0.5, -0.30902, -0.80902, -0.30902, 0.80902, -0.5,
1926 0.30902, -0.80902, 0.5, -0.80902, -0.5, -0.30902, 0.5, -0.30902, -0.80902,
1927 -0.5, 0.30902, 0.80902, -0.30902, 0.80902, -0.5, -0.80902, -0.5, -0.30902,
1928 0, 0, 1, -1, 0, 0, 0, -1, 0,
1929 0.5, -0.30902, 0.80902, 0.30902, -0.80902, -0.5, 0.80902, 0.5, -0.30902,
1930 0.30902, 0.80902, 0.5, 0.80902, -0.5, 0.30902, 0.5, 0.30902, -0.80902,
1931 0.80902, -0.5, 0.30902, -0.5, -0.30902, 0.80902, -0.30902, -0.80902, -0.5,
1932 -0.5, -0.30902, 0.80902, 0.30902, 0.80902, 0.5, -0.80902, 0.5, -0.30902,
1933 0.5, 0.30902, 0.80902, -0.30902, -0.80902, 0.5, 0.80902, -0.5, -0.30902
1934 };
1935
1936 int idx = n % 60;
1937
1938 std::vector<float> matrix(12, 0);
1939 for (int r = 0; r < 3; ++r) {
1940 for (int c = 0; c < 3; ++c) {
1941 matrix[r*4 + c] = matrices[idx*9 + r*3 + c];
1942 }
1943 }
1944
1945 Transform t3d(matrix);
1946 return t3d;
1947}
1948
1949Transform OctahedralSym::get_sym(const int n) const
1950{
1951 // These rotations courtesy of Phil Baldwin
1952 // We have placed the OCT symmetry group with a face along the z-axis
1953 float lvl0 = 0.0f;
1954 float lvl1 = 90.0f;
1955 float lvl2 = 180.0f;
1956
1957 float OCT[72] = {// This is with a face of a cube along z
1958 0,lvl0,0, 0,lvl0,90, 0,lvl0,180, 0,lvl0,270,
1959 0,lvl1,0, 0,lvl1,90, 0,lvl1,180, 0,lvl1,270,
1960 90,lvl1,0, 90,lvl1,90, 90,lvl1,180, 90,lvl1,270,
1961 180,lvl1,0, 180,lvl1,90, 180,lvl1,180, 180,lvl1,270,
1962 270,lvl1,0, 270,lvl1,90, 270,lvl1,180, 270,lvl1,270,
1963 0,lvl2,0, 0,lvl2,90, 0,lvl2,180, 0,lvl2,270
1964 };
1965
1966 int idx = n % 24;
1967// Transform3D ret;
1968// ret.set_rotation((float)OCT[idx * 3 ],(float)OCT[idx * 3 + 1], (float)OCT[idx * 3 + 2] );
1969 Dict d("type","eman");
1970 d["az"] = OCT[idx * 3 ];
1971 d["alt"] = OCT[idx * 3 + 1];
1972 d["phi"] = OCT[idx * 3 + 2];
1973 return Transform(d);
1974
1975}
1976
1977float TetrahedralSym::get_az_alignment_offset() const { return 0.0; }
1978
1979bool TetrahedralSym::is_in_asym_unit(const float& altitude, const float& azimuth, const bool inc_mirror) const
1980{
1981 Dict d = get_delimiters(inc_mirror);
1982 float alt_max = d["alt_max"];
1983 float az_max = d["az_max"];
1984 if ( altitude >= 0 && altitude <= alt_max && azimuth < az_max && azimuth >= 0) {
1985 // convert azimuth to radians
1986 float tmpaz = (float)(EMConsts::deg2rad * azimuth);
1987
1988 float cap_sig = platonic_params["az_max"];
1989 float alt_max = platonic_params["alt_max"];
1990 tmpaz = Util::get_min( tmpaz, cap_sig-tmpaz);
1991
1992 // convert altitude to radians
1993 float tmpalt = (float)(EMConsts::deg2rad * altitude);
1994 if ( platonic_alt_lower_bound(tmpaz, alt_max ) > tmpalt ) {
1995 if ( inc_mirror ) return true;
1996 else {
1997 // you could change the "<" to a ">" here to get the other mirror part of the asym unit
1998 if ( platonic_alt_lower_bound( tmpaz, alt_max/2.0f) < tmpalt ) return false;
1999 else return true;
2000 }
2001 } else return false;
2002 } else return false;
2003}
2004
2005
2006Transform TetrahedralSym::get_sym(const int n) const
2007{
2008 // These rotations courtesy of Phil Baldwin
2009 // It has n=m=3; F=4, E=6=nF/2, V=4=nF/m
2010 float lvl0 = 0.0f; // There is a face along z
2011 float lvl1 = (float)(EMConsts::rad2deg*acos(-1.0f/3.0f)); // There are 3 faces at this angle
2012
2013 float TET[36] = {// This is with the face along z
2014 0,lvl0,0, 0,lvl0,120, 0,lvl0,240,
2015 0,lvl1,60, 0,lvl1,180, 0,lvl1,300,
2016 120,lvl1,60, 120,lvl1,180, 120,lvl1,300,
2017 240,lvl1,60, 240,lvl1,180, 240,lvl1,300
2018 };
2019
2020 int idx = n % 12;
2021
2022 Dict d("type","eman");
2023 d["az"] = TET[idx * 3 ];
2024 d["alt"] = TET[idx * 3 + 1];
2025 d["phi"] = TET[idx * 3 + 2];
2026 return Transform(d);
2027
2028}
2029
2030
2031vector<Vec3f> TetrahedralSym::get_asym_unit_points(bool inc_mirror) const
2032{
2033 vector<Vec3f> ret;
2034
2035 Vec3f b = Vec3f(0,0,1);
2036 ret.push_back(b);
2037 float theta_c_on_two = (float)platonic_params["theta_c_on_two"]; // already in radians
2038 float theta_c = 2*theta_c_on_two;
2039
2040 Vec3f c_on_two = Vec3f(0,-sin(theta_c_on_two),cos(theta_c_on_two));
2041 Vec3f c = Vec3f(0,-sin(theta_c),cos(theta_c));
2042 ret.push_back(c_on_two);
2043 float cap_sig = platonic_params["az_max"];
2044 if ( inc_mirror ) {
2045 Vec3f a = Vec3f(sin(theta_c)*sin(cap_sig),-sin(theta_c)*cos(cap_sig),cos(theta_c));
2046
2047 Vec3f f = a+b+c;
2048 f.normalize();
2049
2050 ret.push_back(f);
2051 }
2052
2053 Vec3f a_on_two = Vec3f(sin(theta_c_on_two)*sin(cap_sig),-sin(theta_c_on_two)*cos(cap_sig),cos(theta_c_on_two));
2054 ret.push_back(a_on_two);
2055
2056
2057 if ( get_az_alignment_offset() != 0 ) {
2058 Dict d("type","eman");
2059 d["az"] = get_az_alignment_offset();
2060 d["phi"] = 0.0f;
2061 d["alt"] = 0.0f;
2062 Transform t(d);
2063 for (vector<Vec3f>::iterator it = ret.begin(); it != ret.end(); ++it )
2064 {
2065 *it = (*it)*t;
2066 }
2067 }
2068
2069 return ret;
2070}
2071
2072
2073
EMAN::CSym
An encapsulation of cyclic 3D symmetry.
Definition: symmetry.h:234
EMAN::CSym::get_delimiters
virtual Dict get_delimiters(const bool inc_mirror=false) const
Get the altitude and phi angle of the c symmetry, which depends on nysm.
Definition: symmetry.cpp:1249
EMAN::CSym::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const
A function to be used when generating orientations over portion of the unit sphere defined by paramet...
Definition: symmetry.cpp:1263
EMAN::CSym::get_asym_unit_triangles
virtual vector< vector< Vec3f > > get_asym_unit_triangles(bool inc_mirror) const
Get triangles that precisely occlude the projection area of the default asymmetric unit.
Definition: symmetry.cpp:1275
EMAN::CSym::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror=false) const
to demarcate the asymmetric unit.
Definition: symmetry.cpp:1370
EMAN::CSym::get_sym
virtual Transform get_sym(const int n) const
Provides access to the complete set of rotational symmetry operations associated with this symmetry.
Definition: symmetry.cpp:1407
EMAN::DSym::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const
A function to be used when generating orientations over portion of the unit sphere defined by paramet...
Definition: symmetry.cpp:1435
EMAN::DSym::get_asym_unit_triangles
virtual vector< vector< Vec3f > > get_asym_unit_triangles(bool inc_mirror) const
Get triangles that precisely occlude the projection area of the default asymmetric unit.
Definition: symmetry.cpp:1472
EMAN::DSym::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror=false) const
Definition: symmetry.cpp:1527
EMAN::DSym::get_delimiters
virtual Dict get_delimiters(const bool inc_mirror=false) const
Get the altitude and phi angle of the d symmetry, which depends on nysm.
Definition: symmetry.cpp:1420
EMAN::DSym::get_sym
virtual Transform get_sym(const int n) const
Provides access to the complete set of rotational symmetry operations associated with this symmetry.
Definition: symmetry.cpp:1452
EMAN::Dict
Dict is a dictionary to store <string, EMObject> pair.
Definition: emobject.h:385
EMAN::Dict::set_default
type set_default(const string &key, type val)
Default setting behavior This can be achieved using a template - d.woolford Jan 2008 (before there wa...
Definition: emobject.h:569
EMAN::EMConsts::rad2deg
static const double rad2deg
Definition: emobject.h:77
EMAN::EMConsts::deg2rad
static const double deg2rad
Definition: emobject.h:78
EMAN::EmanOrientationGenerator::get_az_delta
float get_az_delta(const float &delta, const float &altitude, const int maxcsym) const
Gets the optimum azimuth delta (angular step) for a given altitude, delta and maximum symmetry.
Definition: symmetry.cpp:301
EMAN::EmanOrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
Definition: symmetry.cpp:323
EMAN::EmanOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
generate orientations given some symmetry type
Definition: symmetry.cpp:389
EMAN::EvenOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
Generate even distributed orientations in the asymmetric unit of the symmetry.
Definition: symmetry.cpp:634
EMAN::EvenOrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
Definition: symmetry.cpp:591
EMAN::FactoryBase::copy_relevant_params
Dict copy_relevant_params(const FactoryBase *const that) const
Definition: emobject.h:946
EMAN::FactoryBase::set_params
void set_params(const Dict &new_params)
Set new parameters.
Definition: emobject.h:916
EMAN::FactoryBase::params
Dict params
This is the dictionary the stores the parameters of the object.
Definition: emobject.h:953
EMAN::FactoryBase::get_name
virtual string get_name() const =0
Get the unique name of this class (especially for factory based instantiation access)
EMAN::FactoryBase::get_params
Dict get_params() const
get a copy of the parameters of this class
Definition: emobject.h:911
EMAN::Factory
Factory is used to store objects to create new instances.
Definition: emobject.h:725
EMAN::Factory::get
static T * get(const string &instance_name)
Definition: emobject.h:781
EMAN::HSym::get_delimiters
virtual Dict get_delimiters(const bool inc_mirror=false) const
Get the altitude and phi angle of the d symmetry, which depends on nysm.
Definition: symmetry.cpp:1567
EMAN::HSym::get_asym_unit_triangles
virtual vector< vector< Vec3f > > get_asym_unit_triangles(bool inc_mirror) const
Get triangles that precisely occlude the projection area of the default asymmetric unit.
Definition: symmetry.cpp:1602
EMAN::HSym::get_sym
virtual Transform get_sym(const int n) const
Provides access to the complete set of rotational and translational symmetry operations associated wi...
Definition: symmetry.cpp:1648
EMAN::HSym::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror=false) const
Definition: symmetry.cpp:1608
EMAN::HSym::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const
A function to be used when generating orientations over portion of the unit sphere defined by paramet...
Definition: symmetry.cpp:1585
EMAN::Icosahedral2Sym::get_sym
virtual Transform get_sym(const int n) const
This function provides access to the unique rotational symmetries of an icosahedron.
Definition: symmetry.cpp:1871
EMAN::Icosahedral2Sym::get_az_alignment_offset
virtual float get_az_alignment_offset() const
Get the azimuth alignment offset required to ensure that orientations align correctly with symmetric ...
Definition: symmetry.cpp:1869
EMAN::IcosahedralSym::NAME
static const string NAME
The name of this class - used to access it from factories etc. Should be "icos".
Definition: symmetry.h:910
EMAN::IcosahedralSym::get_sym
virtual Transform get_sym(const int n) const
This function provides access to the unique rotational symmetries of an icosahedron.
Definition: symmetry.cpp:1821
EMAN::IcosahedralSym::get_az_alignment_offset
virtual float get_az_alignment_offset() const
Get the azimuth alignment offset required to ensure that orientations align correctly with symmetric ...
Definition: symmetry.cpp:1819
EMAN::OctahedralSym::get_sym
virtual Transform get_sym(const int n) const
This function provides access to the unique rotational symmetries of an octahedron.
Definition: symmetry.cpp:1949
EMAN::OctahedralSym::NAME
static const string NAME
The name of this class - used to access it from factories etc. Should be "oct".
Definition: symmetry.h:834
EMAN::OptimumOrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
Definition: symmetry.cpp:779
EMAN::OptimumOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
Generate Saff orientations in the asymmetric unit of the symmetry.
Definition: symmetry.cpp:794
EMAN::OptimumOrientationGenerator::optimize_distances
vector< Vec3f > optimize_distances(const vector< Transform > &v) const
Optimize the distances in separating points on the unit sphere, as described by the the rotations in ...
Definition: symmetry.cpp:844
EMAN::OrientationGenerator
An orientation generator is a kind of class that will generate orientations for a given symmetry If o...
Definition: symmetry.h:1000
EMAN::OrientationGenerator::add_orientation
bool add_orientation(vector< Transform > &v, const float &az, const float &alt) const
This functions adds one or more Transform objects to the vector v, depending on the parameters stored...
Definition: symmetry.cpp:273
EMAN::OrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const =0
generate orientations given some symmetry type
EMAN::OrientationGenerator::get_optimal_delta
float get_optimal_delta(const Symmetry3D *const sym, const int &n) const
This function gets the optimal value of the delta (or angular spacing) of the orientations based on a...
Definition: symmetry.cpp:236
EMAN::OrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const =0
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
EMAN::OrientationGenerator::get_az_max
void get_az_max(const Symmetry3D *const sym, const float &altmax, const bool inc_mirror, const float &alt_iterator, const float &h, bool &d_odd_mirror_flag, float &azmax_adjusted) const
Definition: symmetry.cpp:187
EMAN::PlatonicSym::platonic_alt_lower_bound
float platonic_alt_lower_bound(const float &azimuth, const float &alpha) const
Returns the lower bound of the asymmetric unit, as dependent on azimuth, and on alpha - alpha is alt_...
Definition: symmetry.cpp:1737
EMAN::PlatonicSym::get_asym_unit_triangles
virtual vector< vector< Vec3f > > get_asym_unit_triangles(bool inc_mirror) const
Get triangles that precisely occlude the projection area of the default asymmetric unit.
Definition: symmetry.cpp:1750
EMAN::PlatonicSym::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror=false) const
Definition: symmetry.cpp:1777
EMAN::PlatonicSym::init
void init()
Init - Called to initialize platonic_params, should be called in the constructor of all Platonic soli...
Definition: symmetry.cpp:1673
EMAN::PlatonicSym::get_delimiters
virtual Dict get_delimiters(const bool inc_mirror=false) const
Returns the range of altitude and azimuth angles which encompass the asymmetric unit of the Platonic ...
Definition: symmetry.cpp:1694
EMAN::PlatonicSym::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const
A function to be used when generating orientations over portion of the unit sphere defined by paramet...
Definition: symmetry.cpp:1711
EMAN::PlatonicSym::platonic_params
Dict platonic_params
A dictionary that stores important angles, in radians.
Definition: symmetry.h:631
EMAN::RandomOrientationGenerator::NAME
static const string NAME
The name of this class - used to access it from factories etc.
Definition: symmetry.h:1278
EMAN::RandomOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
Generate random orientations in the asymmetric unit of the symmetry.
Definition: symmetry.cpp:546
EMAN::SaffOrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
Definition: symmetry.cpp:681
EMAN::SaffOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
Generate Saff orientations in the asymmetric unit of the symmetry.
Definition: symmetry.cpp:722
EMAN::SingleOrientationGenerator::gen_orientations
virtual vector< Transform > gen_orientations(const Symmetry3D *const sym) const
generate orientations given some symmetry type
Definition: symmetry.cpp:524
EMAN::SingleOrientationGenerator::get_orientations_tally
virtual int get_orientations_tally(const Symmetry3D *const sym, const float &delta) const
This function returns how many orientations will be generated for a given delta (angular spacing) It ...
Definition: symmetry.cpp:519
EMAN::Symmetry3D
Symmetry3D - A base class for 3D Symmetry objects.
Definition: symmetry.h:57
EMAN::Symmetry3D::num_triangles
int num_triangles
This is stores the number of triangles returned by get_asym_unit_triangles(true)
Definition: symmetry.h:212
EMAN::Symmetry3D::au_sym_triangles
vector< vector< Vec3f > > au_sym_triangles
This cache is of size cache_size.
Definition: symmetry.h:214
EMAN::Symmetry3D::cache_au_planes
void cache_au_planes() const
Establish the asymmetric unit planes cache.
Definition: symmetry.cpp:1020
EMAN::Symmetry3D::is_c_sym
virtual bool is_c_sym() const
A function that is used to determine if this is a c symmetry object This function is only virtually o...
Definition: symmetry.h:106
EMAN::Symmetry3D::is_h_sym
virtual bool is_h_sym() const
A function that is used to determine if this is a Helical symmetry object This function is only virtu...
Definition: symmetry.h:100
EMAN::Symmetry3D::get_max_csym
virtual int get_max_csym() const =0
The Symmetry3D object must return the maximum degree of symmetry it exhibits about any one axis.
EMAN::Symmetry3D::is_tet_sym
virtual bool is_tet_sym() const
A function that is used to determine if this is the tetrahedral symmetry object This function is only...
Definition: symmetry.h:118
EMAN::Symmetry3D::get_touching_au_transforms
virtual vector< Transform > get_touching_au_transforms(bool inc_mirror=true) const
Gets a vector of Transform objects that define the set of asymmetric units that touch the default asy...
Definition: symmetry.cpp:1155
EMAN::Symmetry3D::get_nsym
virtual int get_nsym() const =0
The total number of unique symmetry operations that will be return by this object when a calling prog...
EMAN::Symmetry3D::get_asym_unit_triangles
virtual vector< vector< Vec3f > > get_asym_unit_triangles(bool inc_mirror) const =0
Get triangles that precisely occlude the projection area of the default asymmetric unit.
EMAN::Symmetry3D::~Symmetry3D
virtual ~Symmetry3D()
Definition: symmetry.cpp:963
EMAN::Symmetry3D::get_symmetries
static vector< Transform > get_symmetries(const string &symmetry)
Definition: symmetry.cpp:1240
EMAN::Symmetry3D::in_which_asym_unit
virtual int in_which_asym_unit(const Transform &t3d) const
A function that will determine in which asymmetric unit a given orientation resides The asymmetric un...
Definition: symmetry.cpp:1001
EMAN::Symmetry3D::cache_size
int cache_size
Have to remember the cache size.
Definition: symmetry.h:210
EMAN::Symmetry3D::is_platonic_sym
virtual bool is_platonic_sym() const
A function that is used to determine if this is a platonic symmetry object This function is only virt...
Definition: symmetry.h:94
EMAN::Symmetry3D::delete_au_planes
void delete_au_planes()
Clear the asymmetric unit planes cache.
Definition: symmetry.cpp:1059
EMAN::Symmetry3D::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror) const =0
The Symmetry3D object must be capable of returning an ordered list of points on the unit sphere that ...
EMAN::Symmetry3D::point_in_which_asym_unit
virtual int point_in_which_asym_unit(const Vec3f &v) const
A function that will determine in which asymmetric unit a given vector resides The asymmetric unit 'n...
Definition: symmetry.cpp:1072
EMAN::Symmetry3D::get_syms
virtual vector< Transform > get_syms() const
Definition: symmetry.cpp:1224
EMAN::Symmetry3D::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const =0
A function to be used when generating orientations over portion of the unit sphere defined by paramet...
EMAN::Symmetry3D::cached_au_planes
float ** cached_au_planes
The asymmetric unit planes are cached to provide a great speed up the point_in_which_asym_unit functi...
Definition: symmetry.h:207
EMAN::Symmetry3D::get_sym
virtual Transform get_sym(const int n) const =0
Every Symmetry3D object must provide access to the full set of its symmetry operators via this functi...
EMAN::Symmetry3D::ncit
vector< vector< Vec3f > >::iterator ncit
Definition: symmetry.h:60
EMAN::Symmetry3D::get_delimiters
virtual Dict get_delimiters(const bool inc_mirror=false) const =0
Every Symmetry3D object must return a dictionary containing the delimiters that define its asymmetric...
EMAN::Symmetry3D::get_az_alignment_offset
virtual float get_az_alignment_offset() const
This functionality is only relevant to platonic symmetries.
Definition: symmetry.h:86
EMAN::Symmetry3D::is_d_sym
virtual bool is_d_sym() const
A function that is used to determine if this is a d symmetry object This function is only virtually o...
Definition: symmetry.h:112
EMAN::TetrahedralSym::get_asym_unit_points
virtual vector< Vec3f > get_asym_unit_points(bool inc_mirror=false) const
Definition: symmetry.cpp:2031
EMAN::TetrahedralSym::get_az_alignment_offset
virtual float get_az_alignment_offset() const
Get the azimuth alignment offset required to ensure that orientations align correctly with symmetric ...
Definition: symmetry.cpp:1977
EMAN::TetrahedralSym::get_sym
virtual Transform get_sym(const int n) const
This function provides access to the unique rotational symmetries of a tetrahedron.
Definition: symmetry.cpp:2006
EMAN::TetrahedralSym::is_in_asym_unit
virtual bool is_in_asym_unit(const float &altitude, const float &azimuth, const bool inc_mirror) const
In tetrahedral symmetry special consideration must be taken when generating orientations in the asymm...
Definition: symmetry.cpp:1979
EMAN::Transform
A Transform object is a somewhat specialized object designed specifically for EMAN2/Sparx storage of ...
Definition: transform.h:75
EMAN::Transform::ERR_LIMIT
static const float ERR_LIMIT
Definition: transform.h:78
EMAN::Transform::set_trans
void set_trans(const float &x, const float &y, const float &z=0)
Set the post translation component.
Definition: transform.cpp:1036
EMAN::Transform::invert
void invert()
Get the inverse of this transformation matrix.
Definition: transform.cpp:1293
EMAN::Util::str_to_lower
static string str_to_lower(const string &s)
Return a lower case version of the argument string.
Definition: util.cpp:283
EMAN::Util::get_frand
static float get_frand(int low, int high)
Get a float random number between low and high, [low, high)
Definition: util.cpp:725
EMAN::Util::get_min
static int get_min(int f1, int f2)
Get the minimum of 2 numbers.
Definition: util.h:922
EMAN::Util::equation_of_plane
static void equation_of_plane(const Vec3f &p1, const Vec3f &p2, const Vec3f &p3, float *plane)
Determine the equation of a plane that intersects 3 points in 3D space.
Definition: util.cpp:1293
EMAN::Util::get_gauss_rand
static float get_gauss_rand(float mean, float sigma)
Get a Gaussian random number.
Definition: util.cpp:845
EMAN::Vec3::squared_length
Type squared_length() const
Calculate its squared length.
Definition: vec3.h:362
EMAN::Vec3::dot
Type dot(const Vec3< Type2 > &v) const
Calculate the dot product of 'this' vector with a second vector.
Definition: vec3.h:373
EMAN::Vec3::normalize
float normalize()
Normalize the vector and return its length before the normalization.
Definition: vec3.h:332
sqrt
EMData * sqrt() const
return square root of current image
Definition: emdata_core.cpp:1084
InvalidParameterException
#define InvalidParameterException(desc)
Definition: exception.h:361
InvalidValueException
#define InvalidValueException(val, desc)
Definition: exception.h:285
NotExistingObjectException
#define NotExistingObjectException(objname, desc)
Definition: exception.h:130
UnexpectedBehaviorException
#define UnexpectedBehaviorException(desc)
Definition: exception.h:400
EMAN::Symmetry3D::reduce
virtual Transform reduce(const Transform &t3d, int n=0) const
A function that will reduce an orientation, as characterized by Euler anges, into a specific asymmetr...
Definition: symmetry.cpp:974
ENTERFUNC
#define ENTERFUNC
Definition: log.h:48
EXITFUNC
#define EXITFUNC
Definition: log.h:49
EMAN
E2Exception class.
Definition: aligner.h:40
EMAN::Vec3f
Vec3< float > Vec3f
Definition: vec3.h:693
EMAN::dump_orientgens
void dump_orientgens()
Dumps useful information about the OrientationGenerator factory.
Definition: symmetry.cpp:157
EMAN::dump_symmetries_list
map< string, vector< string > > dump_symmetries_list()
dump_symmetries_list, useful for obtaining symmetry information
Definition: symmetry.cpp:70
EMAN::dump_orientgens_list
map< string, vector< string > > dump_orientgens_list()
Can be used to get useful information about the OrientationGenerator factory.
Definition: symmetry.cpp:162
EMAN::dump_symmetries
void dump_symmetries()
dump symmetries, useful for obtaining symmetry information
Definition: symmetry.cpp:65
y
#define y(i, j)
Definition: projector.cpp:1516
x
#define x(i)
Definition: projector.cpp:1517
verify
void verify(const Vec3f &tmp, float *plane, const string &message)
Definition: symmetry.cpp:969
Generated by doxygen 1.9.2