Skip to content

Commit

Permalink
Switch to eigen operations for dot products.
Browse files Browse the repository at this point in the history
  • Loading branch information
morriscb committed Nov 6, 2021
1 parent 23b9f9f commit a2e0de3
Showing 1 changed file with 11 additions and 22 deletions.
33 changes: 11 additions & 22 deletions src/pessimisticPatternMatcherUtils.cc
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,7 @@

#include <cmath>
#include <iostream>
#include "ndarray/eigen.h"
#include "lsst/meas/astrom/pessimisticPatternMatcherUtils.h"

namespace lsst {
Expand All @@ -40,26 +41,15 @@ int check_spoke(double cos_theta_src, double sin_theta_src, ndarray::Array<doubl
for (auto idx = ref_dist_idx_array.begin(); idx != ref_dist_idx_array.end(); idx++) {
/// Compute the delta vector from the pattern center.
unsigned int ref_id = ref_id_array[*idx];
double ref_delta[3];
double proj_ref_delta[3];
ref_delta[0] = reference_array[ref_id][0] - ref_ctr[0];
ref_delta[1] = reference_array[ref_id][1] - ref_ctr[1];
ref_delta[2] = reference_array[ref_id][2] - ref_ctr[2];
ndarray::Array<double, 1, 1> ref_delta = copy(reference_array[ref_id] - ref_ctr);

double ref_dot = ref_delta[0] * ref_ctr[0] + ref_delta[1] * ref_ctr[1] + ref_delta[2] * ref_ctr[2];
proj_ref_delta[0] = ref_delta[0] - ref_dot * ref_ctr[0];
proj_ref_delta[1] = ref_delta[1] - ref_dot * ref_ctr[1];
proj_ref_delta[2] = ref_delta[2] - ref_dot * ref_ctr[2];
double ref_dot = ndarray::asEigenMatrix(ref_delta).dot(ndarray::asEigenMatrix(ref_ctr));
ndarray::Array<double, 1, 1> proj_ref_delta = copy(ref_delta - ref_dot * ref_ctr);
/// Compute the cos between our "center" reference vector and the
/// current reference candidate.
double proj_delta_dist_sq =
(proj_ref_delta[0] * proj_ref_delta[0] + proj_ref_delta[1] * proj_ref_delta[1] +
proj_ref_delta[2] * proj_ref_delta[2]);
double proj_delta_dist_sq = ndarray::asEigenMatrix(proj_ref_delta).dot(ndarray::asEigenMatrix(proj_ref_delta));
double geom_dist_ref = sqrt(proj_ref_ctr_dist_sq * proj_delta_dist_sq);
double cos_theta_ref =
(proj_ref_delta[0] * proj_ref_ctr_delta[0] + proj_ref_delta[1] * proj_ref_ctr_delta[1] +
proj_ref_delta[2] * proj_ref_ctr_delta[2]) /
geom_dist_ref;
double cos_theta_ref = ndarray::asEigenMatrix(proj_ref_delta).dot(ndarray::asEigenMatrix(proj_ref_ctr_delta)) / geom_dist_ref;

/// Make sure we can safely make the comparison in case
/// our "center" and candidate vectors are mostly aligned.
Expand All @@ -85,8 +75,8 @@ int check_spoke(double cos_theta_src, double sin_theta_src, ndarray::Array<doubl
(proj_ref_delta[1] * proj_ref_ctr_delta[2] - proj_ref_delta[2] * proj_ref_ctr_delta[1]) /
geom_dist_ref;
cross_ref[1] =
(proj_ref_delta[2] * proj_ref_ctr_delta[0] - proj_ref_delta[0] * proj_ref_ctr_delta[2]) /
geom_dist_ref;
(proj_ref_delta[2] * proj_ref_ctr_delta[0] - proj_ref_delta[0] * proj_ref_ctr_delta[2]) /
geom_dist_ref;
cross_ref[2] =
(proj_ref_delta[0] * proj_ref_ctr_delta[1] - proj_ref_delta[1] * proj_ref_ctr_delta[0]) /
geom_dist_ref;
Expand All @@ -101,11 +91,10 @@ int check_spoke(double cos_theta_src, double sin_theta_src, ndarray::Array<doubl
sin_comparison = (sin_theta_src - sin_theta_ref) / cos_theta_ref;
}

if (abs(sin_comparison) > src_sin_tol) {
continue;
}
/// Return the correct id of the candidate we found.
return ref_id;
if (abs(sin_comparison) < src_sin_tol) {
return ref_id;
}
}
return -1;
}
Expand Down

0 comments on commit a2e0de3

Please sign in to comment.