79 AutoReg::AutoReg(
Pvl &pvl) {
80 p_template = pvl.
findObject(
"AutoRegistration");
83 p_patternChip.SetSize(3, 3);
84 p_searchChip.SetSize(5, 5);
85 p_fitChip.SetSize(5, 5);
86 p_reducedPatternChip.SetSize(3, 3);
87 p_reducedSearchChip.SetSize(5, 5);
88 p_reducedFitChip.SetSize(5, 5);
89 p_gradientFilterType = None;
91 SetPatternValidPercent(50.0);
92 SetSubsearchValidPercent(50.0);
93 SetPatternZScoreMinimum(1.0);
96 SetSubPixelAccuracy(
true);
97 SetSurfaceModelDistanceTolerance(1.5);
98 SetSurfaceModelWindowSize(5);
100 SetReductionFactor(1);
104 p_totalRegistrations = 0;
105 p_pixelSuccesses = 0;
106 p_subpixelSuccesses = 0;
107 p_patternChipNotEnoughValidDataCount = 0;
108 p_patternZScoreNotMetCount = 0;
109 p_fitChipNoDataCount = 0;
110 p_fitChipToleranceNotMetCount = 0;
111 p_surfaceModelNotEnoughValidDataCount = 0;
112 p_surfaceModelSolutionInvalidCount = 0;
113 p_surfaceModelDistanceInvalidCount = 0;
127 void AutoReg::Init() {
141 for(
int line = 1; line <= p_fitChip.Lines(); line++) {
142 for(
int samp = 1; samp <= p_fitChip.Samples(); samp++) {
149 for(
int line = 1; line <= p_reducedPatternChip.Lines(); line++) {
150 for(
int samp = 1; samp <= p_reducedPatternChip.Samples(); samp++) {
151 p_reducedPatternChip.SetValue(samp, line,
Isis::Null);
157 for(
int line = 1; line <= p_reducedSearchChip.Lines(); line++) {
158 for(
int samp = 1; samp <= p_reducedSearchChip.Samples(); samp++) {
159 p_reducedSearchChip.SetValue(samp, line,
Isis::Null);
166 AutoReg::~AutoReg() {
207 void AutoReg::Parse(
Pvl &pvl) {
211 SetTolerance(algo[
"Tolerance"]);
213 SetChipInterpolator((QString)algo[
"ChipInterpolator"]);
217 SetSubPixelAccuracy((QString)algo[
"SubpixelAccuracy"] ==
"True");
221 SetReductionFactor((
int)algo[
"ReductionFactor"]);
225 SetGradientFilterType((QString)algo[
"Gradient"]);
230 PatternChip()->SetSize((
int)pchip[
"Samples"], (
int)pchip[
"Lines"]);
234 if(pchip.hasKeyword(
"ValidMinimum")) minimum = pchip[
"ValidMinimum"];
235 if(pchip.hasKeyword(
"ValidMaximum")) maximum = pchip[
"ValidMaximum"];
236 PatternChip()->SetValidRange(minimum, maximum);
238 if(pchip.hasKeyword(
"MinimumZScore")) {
239 SetPatternZScoreMinimum((
double)pchip[
"MinimumZScore"]);
241 if(pchip.hasKeyword(
"ValidPercent")) {
242 SetPatternValidPercent((
double)pchip[
"ValidPercent"]);
247 SearchChip()->SetSize((
int)schip[
"Samples"], (
int)schip[
"Lines"]);
251 if(schip.hasKeyword(
"ValidMinimum")) minimum = schip[
"ValidMinimum"];
252 if(schip.hasKeyword(
"ValidMaximum")) maximum = schip[
"ValidMaximum"];
253 SearchChip()->SetValidRange(minimum, maximum);
254 if(schip.hasKeyword(
"SubchipValidPercent")) {
255 SetSubsearchValidPercent((
double)schip[
"SubchipValidPercent"]);
263 SetSurfaceModelDistanceTolerance((
double)smodel[
"DistanceTolerance"]);
267 SetSurfaceModelWindowSize((
int)smodel[
"WindowSize"]);
273 QString msg =
"Improper format for AutoReg PVL [" + pvl.
fileName() +
"]";
286 void AutoReg::SetGradientFilterType(
const QString &gradientFilterType) {
287 if (gradientFilterType ==
"None") {
288 p_gradientFilterType = None;
290 else if (gradientFilterType ==
"Sobel") {
291 p_gradientFilterType = Sobel;
295 "Invalid Gradient type. Cannot use ["
296 + gradientFilterType +
"] to filter chip",
302 QString AutoReg::GradientFilterString()
const {
303 switch (p_gradientFilterType) {
304 case None:
return "None";
305 case Sobel:
return "Sobel";
307 IException::Programmer,
308 "AutoReg only allows Sobel gradient filter or None",
325 void AutoReg::SetSubPixelAccuracy(
bool on) {
326 p_subpixelAccuracy = on;
352 void AutoReg::SetPatternValidPercent(
const double percent) {
353 if((percent <= 0.0) || (percent > 100.0)) {
354 string msg =
"Invalid value for PatternChip ValidPercent ["
356 +
"]. Must be greater than 0.0 and less than or equal to 100.0 (Default is 50.0).";
359 p_patternValidPercent = percent;
380 void AutoReg::SetSubsearchValidPercent(
const double percent) {
381 if((percent <= 0.0) || (percent > 100.0)) {
382 string msg =
"Invalid value for SearchChip SubchipValidPercent ["
384 +
"]. Must be greater than 0.0 and less than or equal to 100.0 (Default is 50.0).";
387 p_subsearchValidPercent = percent;
407 void AutoReg::SetPatternZScoreMinimum(
double minimum) {
409 string msg =
"Invalid value for PatternChip MinimumZScore ["
411 +
"]. Must be greater than 0.0. (Default is 1.0).";
414 p_minimumPatternZScore = minimum;
429 void AutoReg::SetTolerance(
const double tolerance) {
430 p_tolerance = tolerance;
453 void AutoReg::SetChipInterpolator(
const QString &interpolator) {
456 if(interpolator ==
"NearestNeighborType") {
457 itype = Isis::Interpolator::NearestNeighborType;
459 else if(interpolator ==
"BiLinearType") {
460 itype = Isis::Interpolator::BiLinearType;
462 else if(interpolator ==
"CubicConvolutionType") {
463 itype = Isis::Interpolator::CubicConvolutionType;
467 "Invalid Interpolator type. Cannot use ["
468 + interpolator +
"] to load chip",
473 p_patternChip.SetReadInterpolator(itype);
474 p_searchChip.SetReadInterpolator(itype);
475 p_reducedPatternChip.SetReadInterpolator(itype);
476 p_reducedSearchChip.SetReadInterpolator(itype);
493 void AutoReg::SetSurfaceModelWindowSize(
int size) {
494 if(size % 2 != 1 || size < 3) {
495 string msg =
"Invalid value for SurfaceModel WindowSize ["
496 +
IString(size) +
"]. Must be an odd number greater than or equal to 3";
514 void AutoReg::SetSurfaceModelDistanceTolerance(
double distance) {
515 if(distance <= 0.0) {
516 string msg =
"Invalid value for SurfaceModel DistanceTolerance ["
517 +
IString(distance) +
"]. Must greater than 0.0.";
520 p_distanceTolerance = distance;
534 void AutoReg::SetReductionFactor(
int factor) {
536 string msg =
"Invalid value for Algorithm ReductionFactor ["
537 +
IString(factor) +
"]. Must greater than or equal to 1.";
540 p_reduceFactor = factor;
552 Chip AutoReg::Reduce(
Chip &chip,
int reductionFactor) {
554 (int)chip.
Lines() / reductionFactor);
555 if((
int)rChip.Samples() < 1 || (int)rChip.Lines() < 1) {
562 for(
int line = 1; line <= rChip.Lines(); line++) {
563 for(
int samp = 1; samp <= rChip.Samples(); samp++) {
569 for(
int l = 1; l <= rChip.Lines(); l++) {
570 int istartLine = (l - 1) * reductionFactor + 1;
571 int iendLine = istartLine + reductionFactor - 1;
572 for(
int s = 1; s <= rChip.Samples(); s++) {
574 int istartSamp = (s - 1) * reductionFactor + 1;
575 int iendSamp = istartSamp + reductionFactor - 1;
578 for(
int line = istartLine; line < iendLine; line++) {
579 for(
int sample = istartSamp; sample < iendSamp; sample++) {
583 rChip.SetValue(s, l, stats.
Average());
603 int N = p_windowSize / 2 + 1;
605 if(p_searchChip.Samples() < p_patternChip.Samples() + N) {
606 string msg =
"Search chips samples [";
607 msg +=
IString(p_searchChip.Samples()) +
"] must be at ";
608 msg +=
"least [" +
IString(N) +
"] pixels wider than the pattern chip samples [";
609 msg +=
IString(p_patternChip.Samples()) +
"] for successful surface modeling";
613 if(p_searchChip.Lines() < p_patternChip.Lines() + N) {
614 string msg =
"Search chips lines [";
615 msg +=
IString(p_searchChip.Lines()) +
"] must be at ";
616 msg +=
"least [" +
IString(N) +
"] pixels taller than the pattern chip lines [";
617 msg +=
IString(p_patternChip.Lines()) +
"] for successful surface modeling";
622 p_totalRegistrations++;
628 Chip gradientPatternChip(p_patternChip);
629 Chip gradientSearchChip(p_searchChip);
630 ApplyGradientFilter(gradientPatternChip);
631 ApplyGradientFilter(gradientSearchChip);
634 if(!gradientPatternChip.
IsValid(p_patternValidPercent)) {
635 p_patternChipNotEnoughValidDataCount++;
636 p_registrationStatus = PatternChipNotEnoughValidData;
637 return PatternChipNotEnoughValidData;
640 if(!ComputeChipZScore(gradientPatternChip)) {
641 p_patternZScoreNotMetCount++;
642 p_registrationStatus = PatternZScoreNotMet;
643 return PatternZScoreNotMet;
665 int startSamp = (gradientPatternChip.
Samples() - 1) / 2 + 1;
666 int startLine = (gradientPatternChip.
Lines() - 1) / 2 + 1;
667 int endSamp = gradientSearchChip.
Samples() - startSamp + 1;
668 int endLine = gradientSearchChip.
Lines() - startLine + 1;
674 if(gradientPatternChip.
Samples() / p_reduceFactor < 2 || gradientPatternChip.
Lines() / p_reduceFactor < 2) {
675 string msg =
"Reduction factor is too large";
681 int bestSearchSamp = gradientSearchChip.
TackSample();
682 int bestSearchLine = gradientSearchChip.
TackLine();
689 if(p_reduceFactor != 1) {
690 p_reducedPatternChip.SetSize((
int)gradientPatternChip.
Samples() / p_reduceFactor,
691 (int)gradientPatternChip.
Lines() / p_reduceFactor);
696 for(
int line = 1; line <= p_reducedPatternChip.Lines(); line++) {
697 for(
int samp = 1; samp <= p_reducedPatternChip.Samples(); samp++) {
698 p_reducedPatternChip.SetValue(samp, line,
Isis::Null);
702 p_reducedPatternChip =
Reduce(gradientPatternChip, p_reduceFactor);
703 if(!ComputeChipZScore(p_reducedPatternChip)) {
704 p_patternZScoreNotMetCount++;
705 p_registrationStatus = PatternZScoreNotMet;
706 return PatternZScoreNotMet;
709 p_reducedSearchChip =
Reduce(gradientSearchChip, p_reduceFactor);
710 int reducedStartSamp = (p_reducedPatternChip.Samples() - 1) / 2 + 1;
711 int reducedEndSamp = p_reducedSearchChip.Samples() - reducedStartSamp + 1;
712 int reducedStartLine = (p_reducedPatternChip.Lines() - 1) / 2 + 1;
713 int reducedEndLine = p_reducedSearchChip.Lines() - reducedStartLine + 1;
715 Match(p_reducedSearchChip, p_reducedPatternChip, p_reducedFitChip,
716 reducedStartSamp, reducedEndSamp, reducedStartLine, reducedEndLine);
719 p_fitChipNoDataCount++;
720 p_registrationStatus = FitChipNoData;
721 return FitChipNoData;
728 int bs = (p_bestSamp - 1) * p_reduceFactor + ((p_reduceFactor - 1) / 2) + 1;
729 int bl = (p_bestLine - 1) * p_reduceFactor + ((p_reduceFactor - 1) / 2) + 1;
736 int newstartSamp = bs - p_reduceFactor - p_windowSize - 1;
737 int newendSamp = bs + p_reduceFactor + p_windowSize + 1;
738 int newstartLine = bl - p_reduceFactor - p_windowSize - 1;
739 int newendLine = bl + p_reduceFactor + p_windowSize + 1;
741 if(newstartLine < startLine) newstartLine = startLine;
742 if(newendSamp > endSamp) newendSamp = endSamp;
743 if(newstartSamp < startSamp) newstartSamp = startSamp;
744 if(newendLine > endLine) newendLine = endLine;
746 startSamp = newstartSamp;
747 endSamp = newendSamp;
748 startLine = newstartLine;
749 endLine = newendLine;
760 p_registrationStatus = Registration(gradientSearchChip, gradientPatternChip,
761 p_fitChip, startSamp, startLine, endSamp, endLine,
762 bestSearchSamp, bestSearchLine);
765 p_searchChip.SetChipPosition(p_chipSample, p_chipLine);
766 p_cubeSample = gradientSearchChip.
CubeSample();
767 p_cubeLine = gradientSearchChip.
CubeLine();
771 if (p_gradientFilterType != None) {
772 p_gradientSearchChip = gradientSearchChip;
773 p_gradientPatternChip = gradientPatternChip;
776 p_goodnessOfFit = p_bestFit;
779 if (p_registrationStatus == AutoReg::SuccessSubPixel)
780 p_subpixelSuccesses++;
785 return p_registrationStatus;
823 Chip &fChip,
int startSamp,
int startLine,
int endSamp,
int endLine,
824 int bestSamp,
int bestLine) {
827 Match(sChip, pChip, fChip, startSamp, endSamp, startLine, endLine);
832 p_fitChipNoDataCount++;
833 p_registrationStatus = FitChipNoData;
834 return FitChipNoData;
838 if (!CompareFits(p_bestFit, Tolerance())) {
839 p_fitChipToleranceNotMetCount++;
840 p_registrationStatus = FitChipToleranceNotMet;
841 return FitChipToleranceNotMet;
845 if (p_subpixelAccuracy && !IsIdeal(p_bestFit)) {
846 Chip window(p_windowSize, p_windowSize);
847 fChip.
Extract(p_bestSamp, p_bestLine, window);
852 if (!window.
IsValid(100.0 * 2.1 / 3.0)) {
853 p_surfaceModelNotEnoughValidDataCount++;
854 p_registrationStatus = SurfaceModelNotEnoughValidData;
855 return SurfaceModelNotEnoughValidData;
860 bool computedSubPixel = SetSubpixelPosition(window);
861 if (!computedSubPixel)
862 return p_registrationStatus;
866 p_sampMovement = fabs(p_bestSamp - p_chipSample);
867 p_lineMovement = fabs(p_bestLine - p_chipLine);
868 if (p_sampMovement > p_distanceTolerance ||
869 p_lineMovement > p_distanceTolerance) {
871 p_surfaceModelDistanceInvalidCount++;
872 p_registrationStatus = SurfaceModelDistanceInvalid;
873 return SurfaceModelDistanceInvalid;
876 p_registrationStatus = SuccessSubPixel;
879 p_chipSample = p_bestSamp;
880 p_chipLine = p_bestLine;
881 p_registrationStatus = SuccessPixel;
884 return p_registrationStatus;
898 bool AutoReg::ComputeChipZScore(
Chip &chip) {
900 for(
int i = 0; i < chip.
Samples(); i++) {
901 double pixels[chip.
Lines()];
902 for(
int j = 0; j < chip.
Lines(); j++) {
903 pixels[j] = chip.
GetValue(i + 1, j + 1);
915 if (p_zScoreMax < p_minimumPatternZScore && -p_zScoreMin < p_minimumPatternZScore) {
930 void AutoReg::ApplyGradientFilter(
Chip &chip) {
931 if (p_gradientFilterType == None) {
938 if (p_gradientFilterType == Sobel) {
944 "No rule to set sub-chip width for selected Gradient Filter Type.";
955 for (
int line = 1; line <= chip.
Lines(); line++) {
956 for (
int sample = 1; sample <= chip.
Samples(); sample++) {
957 Chip subChip = chip.
Extract(subChipWidth, subChipWidth,
962 Buffer buffer(subChipWidth, subChipWidth, 1, Isis::None);
965 for (
int subChipLine = 1; subChipLine <= subChip.
Lines();
967 for (
int subChipSample = 1; subChipSample <= subChip.
Samples();
969 doubleBuffer[bufferIndex] = subChip.
GetValue(subChipSample,
977 double newPixelValue = 0;
978 if (p_gradientFilterType == Sobel) {
979 SobelGradient(buffer, newPixelValue);
981 filteredChip.SetValue(sample, line, newPixelValue);
986 for (
int line = 1; line <= filteredChip.Lines(); line++) {
987 for (
int sample = 1; sample <= filteredChip.Samples(); sample++) {
988 chip.
SetValue(sample, line, filteredChip.GetValue(sample, line));
1003 void AutoReg::SobelGradient(
Buffer &in,
double &v) {
1004 bool specials =
false;
1005 for(
int i = 0; i < in.
size(); ++i) {
1014 v = abs((in[0] + 2 * in[1] + in[2]) - (in[6] + 2 * in[7] + in[8])) +
1015 abs((in[2] + 2 * in[5] + in[8]) - (in[0] + 2 * in[3] + in[6]));
1034 void AutoReg::Match(
Chip &sChip,
Chip &pChip,
Chip &fChip,
int startSamp,
int endSamp,
int startLine,
int endLine) {
1036 if(startSamp == endSamp && startLine == endLine) {
1037 string msg =
"StartSample [" +
IString(startSamp) +
"] = EndSample ["
1038 +
IString(endSamp) +
"] and StartLine [" +
IString(startLine) +
" = EndLine ["
1046 for(
int line = 1; line <= fChip.
Lines(); line++) {
1047 for(
int samp = 1; samp <= fChip.
Samples(); samp++) {
1055 for(
int line = startLine; line <= endLine; line++) {
1056 for(
int samp = startSamp; samp <= endSamp; samp++) {
1058 sChip.
Extract(samp, line, subsearch);
1061 if(!subsearch.IsValid(p_subsearchValidPercent))
continue;
1064 double fit = MatchAlgorithm(pChip, subsearch);
1069 if((p_bestFit ==
Isis::Null) || CompareFits(fit, p_bestFit)) {
1090 bool AutoReg::SetSubpixelPosition(
Chip &window) {
1093 double samples = window.
Samples();
1094 double lines= window.
Lines();
1096 if (bestDN < window.
GetValue(1, 1)) {
1097 for (
int s=1; s <= samples; s++)
1098 for (
int l=1; l <= lines; l++)
1100 bestDN = 1 / bestDN;
1104 double greatestEdgeDn = 0.0;
1105 for (
int s = 1; s <= samples; s++) {
1106 greatestEdgeDn = max(window.
GetValue(s, 1), greatestEdgeDn);
1107 greatestEdgeDn = max(window.
GetValue(s, lines), greatestEdgeDn);
1109 for (
int l = 2; l <= lines - 1; l++) {
1110 greatestEdgeDn = max(window.
GetValue(1, l), greatestEdgeDn);
1111 greatestEdgeDn = max(window.
GetValue(samples, l), greatestEdgeDn);
1119 double temp = greatestEdgeDn + 0.2 * (bestDN - greatestEdgeDn);
1122 floodFill.setDNRange(temp, 1e100);
1124 Chip selectionChip(window);
1125 floodFill.select(&window, &selectionChip);
1127 double windowSample;
1129 floodFill.centerOfMassWeighted(
1130 &window, &selectionChip, &windowSample, &windowLine);
1132 int offsetS = p_bestSamp - window.
ChipSample();
1133 int offsetL = p_bestLine - window.
ChipLine();
1134 p_chipSample = windowSample + offsetS;
1135 p_chipLine = windowLine + offsetL;
1137 if (p_chipSample != p_chipSample) {
1138 p_surfaceModelSolutionInvalidCount++;
1155 bool AutoReg::CompareFits(
double fit1,
double fit2) {
1156 return(std::fabs(fit1 - IdealFit()) <= std::fabs(fit2 - IdealFit()));
1165 bool AutoReg::IsIdeal(
double fit) {
1166 return(std::fabs(IdealFit() - fit) < 0.00001);
1180 Pvl AutoReg::RegistrationStatistics() {
1182 PvlGroup stats(
"AutoRegStatistics");
1193 PvlGroup grp(
"PatternChipFailures");
1194 grp +=
PvlKeyword(
"PatternNotEnoughValidData",
toString(p_patternChipNotEnoughValidDataCount));
1200 fit +=
PvlKeyword(
"FitChipToleranceNotMet",
toString(p_fitChipToleranceNotMetCount));
1203 PvlGroup model(
"SurfaceModelFailures");
1204 model +=
PvlKeyword(
"SurfaceModelNotEnoughValidData",
toString(p_surfaceModelNotEnoughValidDataCount));
1205 model +=
PvlKeyword(
"SurfaceModelSolutionInvalid",
toString(p_surfaceModelSolutionInvalidCount));
1206 model +=
PvlKeyword(
"SurfaceModelDistanceInvalid",
toString(p_surfaceModelDistanceInvalidCount));
1209 return (AlgorithmStatistics(pvl));
1222 PvlGroup &algo = p_template.findGroup(
"Algorithm", Pvl::Traverse);
1223 reg +=
PvlKeyword(
"Algorithm", algo[
"Name"][0]);
1224 reg +=
PvlKeyword(
"Tolerance", algo[
"Tolerance"][0]);
1226 reg +=
PvlKeyword(
"SubpixelAccuracy", algo[
"SubpixelAccuracy"][0]);
1229 reg +=
PvlKeyword(
"ReductionFactor", algo[
"ReductionFactor"][0]);
1232 reg +=
PvlKeyword(
"Gradient", algo[
"Gradient"][0]);
1235 PvlGroup &pchip = p_template.findGroup(
"PatternChip", Pvl::Traverse);
1236 reg +=
PvlKeyword(
"PatternSamples", pchip[
"Samples"][0]);
1237 reg +=
PvlKeyword(
"PatternLines", pchip[
"Lines"][0]);
1239 reg +=
PvlKeyword(
"PatternMinimum", pchip[
"ValidMinimum"][0]);
1242 reg +=
PvlKeyword(
"PatternMaximum", pchip[
"ValidMaximum"][0]);
1245 reg +=
PvlKeyword(
"MinimumZScore", pchip[
"MinimumZScore"][0]);
1248 SetPatternValidPercent((
double)pchip[
"ValidPercent"]);
1249 reg +=
PvlKeyword(
"ValidPercent", pchip[
"ValidPercent"][0]);
1252 PvlGroup &schip = p_template.findGroup(
"SearchChip", Pvl::Traverse);
1253 reg +=
PvlKeyword(
"SearchSamples", schip[
"Samples"][0]);
1254 reg +=
PvlKeyword(
"SearchLines", schip[
"Lines"][0]);
1256 reg +=
PvlKeyword(
"SearchMinimum", schip[
"ValidMinimum"][0]);
1259 reg +=
PvlKeyword(
"SearchMaximum", schip[
"ValidMaximum"][0]);
1261 if(schip.
hasKeyword(
"SubchipValidPercent")) {
1262 SetSubsearchValidPercent((
double)schip[
"SubchipValidPercent"]);
1263 reg +=
PvlKeyword(
"SubchipValidPercent", schip[
"SubchipValidPercent"][0]);
1266 if(p_template.hasGroup(
"SurfaceModel")) {
1267 PvlGroup &smodel = p_template.findGroup(
"SurfaceModel", Pvl::Traverse);
1269 reg +=
PvlKeyword(
"DistanceTolerance", smodel[
"DistanceTolerance"][0]);
1273 reg +=
PvlKeyword(
"WindowSize", smodel[
"WindowSize"][0]);
1295 reg +=
PvlKeyword(
"Algorithm", AlgorithmName());
1298 SubPixelAccuracy() ?
"True" :
"False");
1300 reg +=
PvlKeyword(
"Gradient", GradientFilterString());
1302 Chip *pattern = PatternChip();
1309 Chip *search = SearchChip();
1315 if (SubPixelAccuracy()) {