@@ -630,6 +630,7 @@ void sys::WireModUtility::ModifyROI(std::vector<float> & roi_data,
630630 double q_orig = 0.0 ;
631631 double q_mod = 0.0 ;
632632 double scale_ratio = 1.0 ;
633+ double sigma_distance = 0.0 ;
633634
634635 // loop over the ticks
635636 for (size_t i_t = 0 ; i_t < roi_data.size (); ++i_t )
@@ -638,6 +639,7 @@ void sys::WireModUtility::ModifyROI(std::vector<float> & roi_data,
638639 q_orig = 0.0 ;
639640 q_mod = 0.0 ;
640641 scale_ratio = 1.0 ;
642+ sigma_distance = 0.0 ;
641643
642644 // loop over the subs
643645 for (auto const & subroi_prop : subROIPropVec)
@@ -647,12 +649,16 @@ void sys::WireModUtility::ModifyROI(std::vector<float> & roi_data,
647649
648650 q_orig += gausFunc (i_t + roi_prop.begin , subroi_prop.center , subroi_prop.sigma , subroi_prop.total_q );
649651 q_mod += gausFunc (i_t + roi_prop.begin , subroi_prop.center , scale_vals.r_sigma * subroi_prop.sigma , scale_vals.r_Q * subroi_prop.total_q );
652+ sigma_distance += ((i_t + roi_prop.begin - subroi_prop.center )**2 / subroi_prop.sigma **2 )*\
653+ gausFunc (i_t + roi_prop.begin , subroi_prop.center , subroi_prop.sigma , subroi_prop.total_q );
650654
651655 if (verbose)
652656 std::cout << " Incrementing q_orig by gausFunc(" << i_t + roi_prop.begin << " , " << subroi_prop.center << " , " << subroi_prop.sigma << " , " << subroi_prop.total_q << " )" << ' \n '
653657 << " Incrementing q_mod by gausFunc(" << i_t + roi_prop.begin << " , " << subroi_prop.center << " , " << scale_vals.r_sigma * subroi_prop.sigma << " , " << scale_vals.r_Q * subroi_prop.total_q << " )" << std::endl;
654658 }
655659
660+ sigma_distance = sigma_distance / q_orig;
661+
656662 // do some sanity checks
657663 if (isnan (q_orig))
658664 {
@@ -664,7 +670,10 @@ void sys::WireModUtility::ModifyROI(std::vector<float> & roi_data,
664670 std::cout << " WARNING: obtained q_mod = NaN... setting scale to 0" << std::endl;
665671 scale_ratio = 0.0 ;
666672 } else if (q_orig < 0.01 ) { // check that this is a sane limit
667- std::cout << " WARNING: obtained q_orig < 0.01 ... setting scale to 1" << std::endl;
673+ if (verbose) std::cout << " WARNING: obtained q_orig < 0.01 ... setting scale to 1" << std::endl;
674+ scale_ratio = 1.0 ;
675+ } else if (sigma_distance > 3 .) {
676+ if (verbose) std::cout << " WARNING: sigma_distance > 3.. setting scale to 1" << std::endl;
668677 scale_ratio = 1.0 ;
669678 } else {
670679 scale_ratio = q_mod / q_orig;
0 commit comments