use combined get_detector_pointing_pixels function instead of generat…

…ing intermediate quaternion vectors
CMB-S4 · Oct 12, 2023 · d77404a · d77404a
1 parent 5581882
commit d77404a
Show file tree

Hide file tree

Showing 8 changed files with 41 additions and 22 deletions.
diff --git a/maps/include/maps/pointing.h b/maps/include/maps/pointing.h
@@ -11,6 +11,10 @@ G3VectorQuat
 get_detector_pointing_quats(double x_offset, double y_offset,
     const G3VectorQuat & trans_quat, MapCoordReference coord_sys);
 
+std::vector<size_t>
+get_detector_pointing_pixels(double x_offset, double y_offset,
+    const G3VectorQuat & trans_quat, G3SkyMapConstPtr skymap);
+
 void get_detector_pointing( 
 	double x_offset, double y_offset,
 	const G3VectorQuat & trans_quat,

diff --git a/maps/src/HitsBinner.cxx b/maps/src/HitsBinner.cxx
@@ -255,11 +255,8 @@ HitsBinner::BinHits(const BolometerProperties &bp, const G3VectorQuat &pointing,
     // Do shared pointers add too much overhead here? Probably not...
 {
 	// Get per-detector pointing timestream
-
-	auto detquats = get_detector_pointing_quats(bp.x_offset, bp.y_offset,
-	    pointing, H->coord_ref);
-
-	auto detpointing = H->QuatsToPixels(detquats);
+	auto detpointing = get_detector_pointing_pixels(bp.x_offset, bp.y_offset,
+	    pointing, H);
 
 	for (size_t i = 0; i < detpointing.size(); i++) {
 		(*H)[detpointing[i]] += 1.0;

diff --git a/maps/src/MapBinner.cxx b/maps/src/MapBinner.cxx
@@ -372,10 +372,8 @@ MapBinner::BinTimestream(const G3Timestream &det, double weight,
 {
 	// Get per-detector pointing timestream
 
-	auto detquats = get_detector_pointing_quats(bp.x_offset, bp.y_offset,
-	    pointing, T->coord_ref);
-
-	auto detpointing = T->QuatsToPixels(detquats);
+	auto detpointing = get_detector_pointing_pixels(bp.x_offset, bp.y_offset,
+	    pointing, T);
 
 	if (Q) {
 		// And polarization coupling

diff --git a/maps/src/MapTODMasker.cxx b/maps/src/MapTODMasker.cxx
@@ -105,9 +105,8 @@ MapTODMasker::Process(G3FramePtr frame, std::deque<G3FramePtr> &out)
 		const BolometerProperties &bp = boloprops_->at(dets[i]);
 
 		// Get per-detector pointing timestream
-		auto detquats = get_detector_pointing_quats(bp.x_offset, bp.y_offset,
-		    *pointing, skymap->coord_ref);
-		auto detpointing = skymap->QuatsToPixels(detquats);
+		auto detpointing = get_detector_pointing_pixels(bp.x_offset, bp.y_offset,
+		    *pointing, skymap);
 
 		det.resize(pointing->size());
 		for (size_t j = 0; j < det.size(); j++)

diff --git a/maps/src/MapTODPointing.cxx b/maps/src/MapTODPointing.cxx
@@ -156,10 +156,8 @@ MapTODPointing::Process(G3FramePtr frame, std::deque<G3FramePtr> &out)
 		const BolometerProperties &bp = boloprops_->at(dets[i]);
 
 		// Get per-detector pointing timestream
-		auto detquats = get_detector_pointing_quats(bp.x_offset, bp.y_offset,
-		    *pointing, stub_->coord_ref);
-
-		auto detpointing = stub_->QuatsToPixels(detquats);
+		auto detpointing = get_detector_pointing_pixels(bp.x_offset, bp.y_offset,
+		    *pointing, stub_);
 
 		for (size_t j = 0; j < sz; j++)
 			det[j] = detpointing[j];

diff --git a/maps/src/SingleDetectorBoresightBinner.cxx b/maps/src/SingleDetectorBoresightBinner.cxx
@@ -184,9 +184,7 @@ SingleDetectorBoresightBinner::Process(G3FramePtr frame,
 	g3_assert(timestreams->NSamples() == pointing->size());
 
 	// Conjugate pointing rotation with boresight vector
-	auto pointing_vec = get_detector_pointing_quats(0, 0,
-	    *pointing, template_->coord_ref);
-	auto pixels = template_->QuatsToPixels(pointing_vec);
+	auto pixels = get_detector_pointing_pixels(0, 0, *pointing, template_);
 
 	for (size_t i = 0; i < pixels.size(); i++)
 		(*map_weights_->TT)[pixels[i]] += 1;

diff --git a/maps/src/SingleDetectorMapBinner.cxx b/maps/src/SingleDetectorMapBinner.cxx
@@ -160,10 +160,9 @@ SingleDetectorMapBinner::Process(G3FramePtr frame,
 		if (bp == boloprops_->end())
 			log_fatal("Missing bolometer properties for %s",
 			    det.c_str());
-		auto detquats = get_detector_pointing_quats(
+		auto pixels = get_detector_pointing_pixels(
 		    bp->second.x_offset, bp->second.y_offset,
-		    *pointing, m->coord_ref);
-		auto pixels = m->QuatsToPixels(detquats);
+		    *pointing, m);
 
 		g3_assert(ts->size() == pixels.size());
 		for (size_t j = 0; j < ts->size(); j++) {

diff --git a/maps/src/pointing.cxx b/maps/src/pointing.cxx
@@ -383,6 +383,32 @@ get_detector_pointing_quats(double x_offset, double y_offset,
 	return det_quats;
 }
 
+std::vector<size_t>
+get_detector_pointing_pixels(double x_offset, double y_offset,
+    const G3VectorQuat &trans_quat, G3SkyMapConstPtr skymap)
+{
+	quat q_off = offsets_to_quat(x_offset, y_offset);
+	size_t nsamp = trans_quat.size();
+	std::vector<size_t> pixels(nsamp, (size_t) -1);
+	quat q;
+
+	if (skymap->coord_ref == Local) {
+		for (size_t i = 0; i < nsamp; i++) {
+			q = trans_quat[i] * q_off * ~trans_quat[i];
+			q = quat(q.R_component_1(), q.R_component_2(),
+			    q.R_component_3(), -q.R_component_4());
+			pixels[i] = skymap->QuatToPixel(q);
+		}
+	} else {
+		for (size_t i = 0; i < nsamp; i++) {
+			q = trans_quat[i] * q_off * ~trans_quat[i];
+			pixels[i] = skymap->QuatToPixel(q);
+		}
+	}
+
+	return pixels;
+}
+
 void
 get_detector_pointing(double x_offset, double y_offset,
     const G3VectorQuat &trans_quat, MapCoordReference coord_sys,