EMAN2
Public Types | Public Member Functions | Static Public Member Functions | Static Private Attributes | List of all members
EMAN::BoxingTools Class Reference

BoxingTools is class for encapsulating common boxing operations that may become expensive if they are implemented in python. More...

#include <boxingtools.h>

Public Types

enum  CmpMode { SWARM_DIFFERENCE , SWARM_RATIO , SWARM_AVERAGE_RATIO }
 

Public Member Functions

 BoxingTools ()
 
 ~BoxingTools ()
 

Static Public Member Functions

static vector< float > get_min_delta_profile (const EMData *const image, int x, int y, int radius)
 Gets a pixel minimum delta radial profile about some pixel focal point. More...
 
static bool is_local_maximum (const EMData *const image, int x, int y, int radius, EMData *const exclusion_map)
 Determines if a given pixel is the maximum in local radial neighborhood Useful for automated correlation-based boxing. More...
 
static vector< IntPointauto_correlation_pick (const EMData *const image, float threshold, int radius, const vector< float > &profile, EMData *const exclusion, const int cradius, int mode=1)
 
static bool hi_brid (const EMData *const image, int x, int y, int radius, EMData *const exclusion_map, vector< float > &profile)
 
static void set_radial_non_zero (EMData *const exclusion, int x, int y, int radius)
 
static IntPoint find_radial_max (const EMData *const map, int x, int y, int radius)
 
static map< unsigned int, unsigned int > classify (const vector< vector< float > > &data, const unsigned int &classes=4)
 
static Vec3f get_color (const unsigned int index)
 
static void set_region (EMData *const image, const EMData *const mask, const int x, const int y, const float &val)
 
static void set_mode (const CmpMode m)
 

Static Private Attributes

static vector< Vec3fcolors = vector<Vec3f>()
 
static CmpMode mode = SWARM_AVERAGE_RATIO
 

Detailed Description

BoxingTools is class for encapsulating common boxing operations that may become expensive if they are implemented in python.

Author
David Woolford
Date
April 2008

Definition at line 56 of file boxingtools.h.

Member Enumeration Documentation

◆ CmpMode

Enumerator
SWARM_DIFFERENCE 
SWARM_RATIO 
SWARM_AVERAGE_RATIO 

Definition at line 102 of file boxingtools.h.

Constructor & Destructor Documentation

◆ BoxingTools()

EMAN::BoxingTools::BoxingTools ( )
inline

Definition at line 59 of file boxingtools.h.

59{}

◆ ~BoxingTools()

EMAN::BoxingTools::~BoxingTools ( )
inline

Definition at line 60 of file boxingtools.h.

60{}

Member Function Documentation

◆ auto_correlation_pick()

vector< IntPoint > BoxingTools::auto_correlation_pick ( const EMData *const  image,
float  threshold,
int  radius,
const vector< float > &  profile,
EMData *const  exclusion,
const int  cradius,
int  mode = 1 
)
static

Definition at line 608 of file boxingtools.cpp.

609{
610 if (mode < 0 || mode > 2 ) {
611 throw InvalidValueException(mode,"Error, the mode can only be 0,1, or 2.");
612 }
613
614 if ( radius < 0) {
615 throw InvalidValueException(radius,"Radius must be greater than 1");
616 }
617
618 if ( cradius < 0) {
619 throw InvalidValueException(cradius,"CRadius must be greater than 1");
620 }
621
622
623 int nx = image->get_xsize();
624 int ny = image->get_ysize();
625
626 vector<IntPoint> solution;
627
628 int r = radius+1;
629
630 for(int j = r; j < ny-r;++j) {
631 for(int k = r; k < nx-r;++k) {
632
633 if (exclusion->get_value_at(k,j) != 0 ) continue;
634
635 if (image->get_value_at(k,j) < threshold) continue;
636 if ( mode == 0 ) {
637 solution.push_back(IntPoint(k,j));
638 set_radial_non_zero(exclusion,k,j,radius);
639 continue;
640 }
641
642 vector<float> p(r,0);
643
644 if (hi_brid(image,k,j,r,exclusion,p)) {
645 if ( mode == 1 ) {
646 if (p[cradius] >= profile[cradius]) {
647 solution.push_back(IntPoint(k,j));
648 set_radial_non_zero(exclusion,k,j,radius);
649 continue;
650 }
651 }
652 else /* mode == 2 */{
653 bool bad = false;
654 for (int ii = 0; ii <= cradius; ++ii) {
655 if (p[ii] < profile[ii]) {
656 bad = true;
657 break;
658 }
659 }
660 if (bad) continue;
661 solution.push_back(IntPoint(k,j));
662 set_radial_non_zero(exclusion,k,j,radius);
663 }
664
665
666 }
667 }
668 }
669 return solution;
670}
static void set_radial_non_zero(EMData *const exclusion, int x, int y, int radius)
static bool hi_brid(const EMData *const image, int x, int y, int radius, EMData *const exclusion_map, vector< float > &profile)
static CmpMode mode
Definition: boxingtools.h:113
IntPoint defines an integer-coordinate point in a 1D/2D/3D space.
Definition: geometry.h:192
#define InvalidValueException(val, desc)
Definition: exception.h:285

References hi_brid(), InvalidValueException, mode, and set_radial_non_zero().

◆ classify()

map< unsigned int, unsigned int > BoxingTools::classify ( const vector< vector< float > > &  data,
const unsigned int &  classes = 4 
)
static

Definition at line 850 of file boxingtools.cpp.

851{
852 BoxSVDClassifier classifier(data, classes);
853 map< unsigned int, unsigned int> mapping = classifier.go();
854
856
857 return mapping;
858
859}
static map< unsigned int, unsigned int > colorMappingByClassSize(const map< unsigned int, unsigned int > &grouping)

References EMAN::BoxSVDClassifier::colorMappingByClassSize(), and EMAN::BoxSVDClassifier::go().

◆ find_radial_max()

IntPoint BoxingTools::find_radial_max ( const EMData *const  map,
int  x,
int  y,
int  radius 
)
static

Definition at line 779 of file boxingtools.cpp.

780{
781 float currentmax = map->get_value_at(x,y);
782
783 IntPoint soln(x,y);
784
785 int radius_squared = radius*radius;
786 for(int k = -radius; k <= radius; ++k) {
787 for(int j = -radius; j <= radius; ++j) {
788 // Translate coordinates
789 int xx = x+j;
790 int yy = y+k;
791
792 // Protect against accessing pixels out of bounds
793 if ( xx >= map->get_xsize() || xx < 0 ) continue;
794 if ( yy >= map->get_ysize() || yy < 0 ) continue;
795
796 // Protect against vector accesses beyond the boundary
797 int square_length = k*k + j*j;
798 if (square_length > radius_squared ) continue;
799
800 float val = map->get_value_at(xx,yy);
801
802 if (val > currentmax) {
803 currentmax = val;
804 soln[0] = xx;
805 soln[1] = yy;
806 }
807 }
808 }
809
810 return soln;
811}
#define y(i, j)
Definition: projector.cpp:1516
#define x(i)
Definition: projector.cpp:1517

References x, and y.

◆ get_color()

Vec3f BoxingTools::get_color ( const unsigned int  index)
static

Definition at line 862 of file boxingtools.cpp.

863{
864 if ( colors.size() == 0 ) {
865 colors.push_back(Vec3f(0,0,0));
866 colors.push_back(Vec3f(0,0,1));
867 colors.push_back(Vec3f(0,1,0));
868 colors.push_back(Vec3f(1,0,0));
869 colors.push_back(Vec3f(1,1,0));
870 colors.push_back(Vec3f(1,0,1));
871 colors.push_back(Vec3f(0,1,1));
872 colors.push_back(Vec3f(1,1,1));
873 }
874 if ( index >= colors.size() )
875 {
876 while ( colors.size() <= index )
877 {
878 bool found = false;
879 while ( !found )
880 {
881 unsigned int random_idx = rand() % colors.size();
882 unsigned int random_idx2 = rand() % colors.size();
883 while ( random_idx2 == random_idx )
884 {
885 random_idx2 = rand() % colors.size();
886 }
887
888 Vec3f result = (colors[random_idx] + colors[random_idx2])/2.0;
889 if ( find( colors.begin(), colors.end(), result ) == colors.end() )
890 {
891 colors.push_back( result );
892 found = true;
893 }
894 }
895 }
896 }
897 return colors[index];
898}
static vector< Vec3f > colors
Definition: boxingtools.h:112
Vec3< float > Vec3f
Definition: vec3.h:693

References colors.

◆ get_min_delta_profile()

vector< float > BoxingTools::get_min_delta_profile ( const EMData *const  image,
int  x,
int  y,
int  radius 
)
static

Gets a pixel minimum delta radial profile about some pixel focal point.

Useful for automated correlation-based boxing.

Parameters
imagethe image containing the interesting pixels values (typically a correlation image)
xthe x coordinate of the pixel focal point
ythe y coordinate of the pixel focal point
radiusthe constraining radius of the minimum delta profile
Author
David Woolford
Date
April 2008

Definition at line 498 of file boxingtools.cpp.

499{
500 float peakval = image->get_value_at(x,y);
501
502 vector<float> profile(radius,0); // this sets the vectors size to radius, and the values to 0
503 int radius_squared = radius*radius;
504
505 static vector<float> squared_numbers;
506 if ( (unsigned int)(radius+1) > squared_numbers.size() ) {
507 for(int i = squared_numbers.size(); i <= radius; ++i) {
508 squared_numbers.push_back((float)(i*i));
509 }
510 }
511
512 vector<unsigned int> tally;
513 if (mode == SWARM_AVERAGE_RATIO) tally.resize(profile.size(),0);
514
515 for(int k = -radius; k <= radius; ++k) {
516 for(int j = -radius; j <= radius; ++j) {
517 // Translate coordinates
518 int xx = x+j;
519 int yy = y+k;
520
521 // Protect against accessing pixels out of bounds
522 if ( xx >= image->get_xsize() || xx < 0 ) continue;
523 if ( yy >= image->get_ysize() || yy < 0 ) continue;
524
525 // We don't need to pay attention to the origin
526 if ( xx == x && yy == y) continue;
527
528 // Protect against vector accesses beyond the boundary
529 int square_length = k*k + j*j;
530 if (square_length > radius_squared ) continue;
531
532 // The idx is the radius, rounded down. This creates a certain type of pattern that
533 // can only really be explained visually...
534 int idx = -1;
535 // This little loop avoids the use of sqrtf
536 for(int i = 1; i < radius; ++i) {
537 if ( square_length >= squared_numbers[i] && square_length <= squared_numbers[i+1] ) {
538 idx = i;
539 }
540 }
541// int idx = (int) sqrtf(k*k + j*j);
542 // decrement the idx, because the origin information is redundant
543 idx -= 1;
544
545 if ( mode == SWARM_DIFFERENCE ) {
546 // Finally, get the drop relative to the origin
547 float val = peakval - image->get_value_at(xx,yy);
548
549 // Store it if the drop is smaller than the current value (or if there is no value)
550 if ( profile[idx] > val || profile[idx] == 0 ) profile[idx] = val;
551 }
552 else if (mode == SWARM_RATIO) {
553 // Finally, get the drop relative to the origin
554 float val = (peakval - image->get_value_at(xx,yy) ) / peakval;
555
556 // Store it if the drop is smaller than the current value (or if there is no value)
557 if ( profile[idx] > val || profile[idx] == 0 ) profile[idx] = val;
558 }
559 else if (mode == SWARM_AVERAGE_RATIO) {
560 profile[idx] += image->get_value_at(xx,yy);
561 tally[idx]++;
562 }
563
564 }
565 }
566
567 if (mode == SWARM_AVERAGE_RATIO) {
568 for(unsigned int i = 0; i < profile.size(); ++i) {
569 if (tally[i] != 0) {
570 profile[i] /= static_cast<float>(tally[i]);
571 profile[i] = (peakval - profile[i] ) / peakval;
572 }
573 }
574 }
575
576 return profile;
577}

References mode, SWARM_AVERAGE_RATIO, SWARM_DIFFERENCE, SWARM_RATIO, x, and y.

◆ hi_brid()

bool BoxingTools::hi_brid ( const EMData *const  image,
int  x,
int  y,
int  radius,
EMData *const  exclusion_map,
vector< float > &  profile 
)
static

Definition at line 673 of file boxingtools.cpp.

674{
675
676 float peakval = image->get_value_at(x,y);
677
678 int radius_squared = radius*radius;
679
680 static vector<float> squared_numbers;
681 if ( (unsigned int)(radius+1) > squared_numbers.size() ) {
682 for(int i = squared_numbers.size(); i <= radius; ++i) {
683 squared_numbers.push_back((float)(i*i));
684 }
685 }
686
687 vector<unsigned int> tally;
688 if (mode == SWARM_AVERAGE_RATIO) tally.resize(profile.size(),0);
689
690 for(int k = -radius; k <= radius; ++k) {
691 for(int j = -radius; j <= radius; ++j) {
692 // Translate coordinates
693 int xx = x+j;
694 int yy = y+k;
695
696 // Protect against accessing pixels out of bounds
697 if ( xx >= image->get_xsize() || xx < 0 ) continue;
698 if ( yy >= image->get_ysize() || yy < 0 ) continue;
699
700 // We don't need to pay attention to the origin
701 if ( xx == x && yy == y) continue;
702
703 // Protect against vector accesses beyond the boundary
704 int square_length = k*k + j*j;
705 if (square_length > radius_squared ) continue;
706
707 // It's not a local maximum!
708 if ( image->get_value_at(xx,yy) > peakval) return false;
709
710 // The idx is the radius, rounded down. This creates a certain type of pattern that
711 // can only really be explained visually...
712 int idx = -1;
713 // This little loop avoids the use of sqrtf
714 for(int i = 1; i < radius; ++i) {
715 if ( square_length >= squared_numbers[i] && square_length <= squared_numbers[i+1] ) {
716 idx = i;
717 }
718 }
719// int idx = (int) sqrtf(k*k + j*j);
720 // decrement the idx, because the origin information is redundant
721 idx -= 1;
722
723 if (mode == SWARM_DIFFERENCE) {
724 // Finally, get the drop relative to the origin
725 float val = peakval - image->get_value_at(xx,yy);
726
727 // Store it if the drop is smaller than the current value (or if there is no value)
728 if ( profile[idx] > val || profile[idx] == 0 ) profile[idx] = val;
729 }
730 else if (mode == SWARM_RATIO) {
731 // Finally, get the drop relative to the origin
732 float val = (peakval - image->get_value_at(xx,yy) ) / peakval;
733
734 // Store it if the drop is smaller than the current value (or if there is no value)
735 if ( profile[idx] > val || profile[idx] == 0 ) profile[idx] = val;
736 }
737 else if (mode == SWARM_AVERAGE_RATIO) {
738 profile[idx] += image->get_value_at(xx,yy);
739 tally[idx]++;
740 }
741
742 }
743 }
744
745 if (mode == SWARM_AVERAGE_RATIO) {
746 for(unsigned int i = 0; i < profile.size(); ++i) {
747 if (tally[i] != 0) {
748 profile[i] /= static_cast<float>(tally[i]);
749 profile[i] = (peakval - profile[i] ) / peakval;
750 }
751 }
752 }
753
754 set_radial_non_zero(exclusion_map,x,y,radius);
755
756 return true;
757}

References mode, set_radial_non_zero(), SWARM_AVERAGE_RATIO, SWARM_DIFFERENCE, SWARM_RATIO, x, and y.

Referenced by auto_correlation_pick().

◆ is_local_maximum()

bool BoxingTools::is_local_maximum ( const EMData *const  image,
int  x,
int  y,
int  radius,
EMData *const  exclusion_map 
)
static

Determines if a given pixel is the maximum in local radial neighborhood Useful for automated correlation-based boxing.

Parameters
imagethe image containing the interesting pixels values (typically a correlation image)
xthe x coordinate of the candidate pixel maximum
ythe y coordinate of the candidate pixel maximum
radiusthe constraining radius of the local neighborhood
exclusion_mapan EMData object with the same dimensions as the main image. If a local maximum is found, the pixels in the radial neighborhood that was queried are set to 0. This can be exploited by the calling function to minimize queries.
Author
David Woolford
Date
April 2008

Definition at line 579 of file boxingtools.cpp.

580{
581 float peakval = image->get_value_at(x,y);
582 int radius_squared = radius*radius;
583 for(int k = -radius; k <= radius; ++k) {
584 for(int j = -radius; j <= radius; ++j) {
585 // Translate coordinates
586 int xx = x+j;
587 int yy = y+k;
588
589// // Protect against accessing pixel out of bounds
590 if ( xx >= image->get_xsize() || xx < 0 ) continue;
591 if ( yy >= image->get_ysize() || yy < 0 ) continue;
592
593 // We don't need to pay attention to the origin
594 if ( xx == x && yy == y) continue;
595
596 if ((k*k+j*j)>radius_squared) continue;
597
598 if ( image->get_value_at(xx,yy) > peakval) return false;
599 }
600 }
601
602 set_radial_non_zero(exclusion_map,x,y,radius);
603
604 return true;
605
606}

References set_radial_non_zero(), x, and y.

◆ set_mode()

static void EMAN::BoxingTools::set_mode ( const CmpMode  m)
inlinestatic

Definition at line 108 of file boxingtools.h.

108{ mode = m; }

References mode.

◆ set_radial_non_zero()

void BoxingTools::set_radial_non_zero ( EMData *const  exclusion,
int  x,
int  y,
int  radius 
)
static

Definition at line 760 of file boxingtools.cpp.

761{
762 int radius_squared = radius*radius;
763 for(int k = -radius; k <= radius; ++k) {
764 for(int j = -radius; j <= radius; ++j) {
765 // Translate coordinates
766 int xx = x+j;
767 int yy = y+k;
768
769 if ((k*k+j*j)>radius_squared) continue;
770 // Protect against accessing pixel out of bounds
771 if ( xx >= exclusion->get_xsize() || xx < 0 ) continue;
772 if ( yy >= exclusion->get_ysize() || yy < 0 ) continue;
773
774 exclusion->set_value_at(xx,yy,1);
775 }
776 }
777}

References x, and y.

Referenced by auto_correlation_pick(), hi_brid(), and is_local_maximum().

◆ set_region()

void BoxingTools::set_region ( EMData *const  image,
const EMData *const  mask,
const int  x,
const int  y,
const float &  val 
)
static

Definition at line 814 of file boxingtools.cpp.

814 {
815
816 // Works only in 2D
817 int inx = image->get_xsize();
818 int iny = image->get_ysize();
819 int mnx = mask->get_xsize();
820 int mny = mask->get_ysize();
821
822 int startx = x-mnx/2;
823 int endx =startx + mnx;
824 int xoffset = 0;
825 if (startx < 0) {
826 xoffset = abs(startx);
827 startx = 0;
828 }
829 if (endx > inx) endx = inx;
830
831 int starty = y-mny/2;
832 int endy =starty + mny;
833 int yoffset = 0;
834 if (starty < 0) {
835 yoffset = abs(starty);
836 starty = 0;
837 }
838 if (endy > iny) endy = iny;
839
840
841 for (int j = starty; j < endy; ++j ) {
842 for (int i = startx; i < endx; ++i) {
843 if (mask->get_value_at(xoffset+i-startx,yoffset+j-starty) != 0 ) {
844 image->set_value_at(i,j,val);
845 }
846 }
847 }
848}

References x, and y.

Member Data Documentation

◆ colors

vector< Vec3f > BoxingTools::colors = vector<Vec3f>()
staticprivate

Definition at line 112 of file boxingtools.h.

Referenced by get_color().

◆ mode

BoxingTools::CmpMode BoxingTools::mode = SWARM_AVERAGE_RATIO
staticprivate

Definition at line 113 of file boxingtools.h.

Referenced by auto_correlation_pick(), get_min_delta_profile(), hi_brid(), and set_mode().


The documentation for this class was generated from the following files: