Fixed problem with apollopanstitcher test not passing on OSX. Fixes #4948 (6f78d8cc) · Commits · aflab / astrogeology / ISIS3

isis/src/apollo/apps/apollopanstitcher/apollopanstitcher.cpp

+270 −227

Original line number	Diff line number	Diff line
		@@ -43,6 +43,9 @@ Program uses the fiducial marks along the edges of Apollo Panoramic Images to st
		* @internal
		* @History 2011-11-02 Orrin Thomas - Original Version.
		*
		* @History 2017-06-29 Christopher Combs - Moved initialization of trans[7] to before
		* the loop that uses it. Fixes #4948.
		*
		*/


		@@ -128,8 +131,12 @@ void IsisMain() {
		averageLines = AVERl *5.0/resolution; //scaled average diference between the top
		// and bottom fiducials

		if (15.0/resolution < 1.5) play=1.5;
		else play = 15.0/resolution;
		if (15.0 / resolution < 1.5) {
		play = 1.5;
		}
		else {
		play = 15.0 / resolution;
		}

		//copy the patternS chip (the entire ApolloPanFiducialMark.cub)
		FileName fiducialFileName("$apollo15/calibration/ApolloPanFiducialMark.cub");
		@@ -169,10 +176,12 @@ void IsisMain() {
		//Centroid centroid;
		CentroidApolloPan centroid(resolution);
		Chip inputChip,selectionChip;
		if( panC[0]->pixelType() == 1) //UnsignedByte
		if (panC[0]->pixelType() == 1) { //UnsignedByte
		centroid.setDNRange(12, 1e99); //8 bit bright target
		else
		}
		else {
		centroid.setDNRange(3500, 1e99); //16 bit bright target
		}
		inputChip.SetSize(int(ceil(200 * 5.0 / resolution)), int(ceil(200 * 5.0 / resolution)));

		//find conjugate finducials and calcualte transformations
		@@ -260,8 +269,9 @@ void IsisMain() {

		//look for the top fiducial
		//is this is the first time through the loop?
		if (s == scanFid[0][0])
		if (s == scanFid[0][0]) {
		continue; //then the top fiducial was already found
		}
		searchS.TackCube(s, l);
		searchS.Load(*panC[i], 0, scale);
		*arS->SearchChip() = searchS;
		@@ -371,16 +381,15 @@ void IsisMain() {

		//clear out previous data
		solV.clear();
		//matrix names loosly following the naming convention of Uotilia see large comment above
		//matrix names loosely following the naming convention of Uotilia see large comment above
		double wdot[2], adot[2][3], ata[6], atf[3];

		int l,m,n,o,ni; //additonal indeces for matrix multipilcations, etc
		int l, m, n, o, ni; //additonal indices for matrix multipilcations, etc

		for (j = 0; j < int(conFid.size()) - 1; j++) {
		for (k = j + 1; k < int(conFid.size()); k++) {
		//Starting by finding the two point solution: using the following close form algorithim
		sol.theta = atan2( conFid[j][3] - conFid[k][3],
		conFid[j][2] - conFid[k][2] ) -
		sol.theta = atan2(conFid[j][3] - conFid[k][3], conFid[j][2] - conFid[k][2]) -
		atan2(conFid[j][1] - conFid[k][1], conFid[j][0] - conFid[k][0]);
		sol.dx = conFid[j][2] - conFid[j][0] * cos(sol.theta) + conFid[j][1] * sin(sol.theta);
		sol.dy = conFid[j][3] - conFid[j][0] * sin(sol.theta) - conFid[j][1] * cos(sol.theta);
		@@ -444,8 +453,9 @@ void IsisMain() {
		//mdot[1] = 0.0;
		//mdot[2] = 0.5;

		if (sqrt(0.5(wdot[0]wdot[0] + wdot[1]*wdot[1])) > 3.0)
		if (sqrt(0.5 * (wdot[0] * wdot[0] + wdot[1] * wdot[1])) > 3.0) {
		continue; //if the R^2 residual is greater than 3.0 pixels go on to the next point
		}

		//partials wrt unknowns
		adot[0][0] = -conFid[m][0]st - conFid[m][1]ct; //wrt to theta
		@@ -460,10 +470,13 @@ void IsisMain() {
		//build the normal equations
		//add transpose(adot)*adot to ata
		//note: becuase m is constant, and every weight is equal, it is irrelevant and ignored
		for (l=0;l<3;l++)
		for (n=0;n<=l;n++)
		for (o=0;o<2;o++)
		for (l = 0; l < 3; l++) {
		for (n = 0; n <= l; n++) {
		for (o = 0; o < 2; o++) {
		ata[isymp(l,n)] += adot[o][l] * adot[o][n];
		}
		}
		}

		//add transpose(adot)*wdot to atf
		//note: becuase m is constant, and every weight is equal, it is irrelevant and ignored
		@@ -474,8 +487,9 @@ void IsisMain() {

		//solve normal equations
		//calculation is done in place as delta is returned in atf
		if (choleski_solve(ata,atf,3,2) != 1)
		if (choleski_solve(ata, atf, 3, 2) != 1) {
		continue; //if the solution fails move on to the next two point solution
		}

		//check for nonsense numbers
		if (atf[2] != atf[2]) continue;
		@@ -489,8 +503,8 @@ void IsisMain() {
		if (fabs(atf[0]) < 1e-10 && fabs(atf[0]) < 1e-5 && fabs(atf[0]) < 1e-5) {
		// solution converged
		// find residual stats
		//Note residuals are caluclated for all points-including any that might have been
		// excluded from thes solution above
		// Note: residuals are caluclated for all points-including any that might have been
		// excluded from the solution above
		for (m = 0, sol.maxR = 0.0, sol.averR = 0.0, sol.limit[0] = 0.0, sol.limit[1] = 0.0;
		m<int(conFid.size());
		m++) { // for each point
		@@ -498,7 +512,9 @@ void IsisMain() {
		wdot[1] = -( conFid[m][0]st + conFid[m][1]ct + sol.dy - conFid[m][3]);
		temp = sqrt(wdot[0]wdot[0] + wdot[1]wdot[1]);

		if (temp > sol.maxR) sol.maxR = temp;
		if (temp > sol.maxR) {
		sol.maxR = temp;
		}

		sol.averR += temp;

		@@ -509,10 +525,9 @@ void IsisMain() {

		// save the solution
		solV.push_back(sol);
		break; //break out of interation counting loop
		break; // break out of iteration counting loop
		}
		} // iteration counting for loop

		}
		} // end of two point solutions loops
		if (solV.size() < 1) {
		@@ -526,27 +541,39 @@ void IsisMain() {
		// first lets thin by maximum residual
		temp = solV[0].maxR; // find the lowest max residual
		for (j = 1; j < int(solV.size()); j++) {
		if (solV[j].maxR < temp) temp = solV[j].maxR;
		if (solV[j].maxR < temp) {
		temp = solV[j].maxR;
		}
		}
		// toss out any solution that has a maximum residual higher than ceil(lowest max residual)
		temp = ceil(temp);
		for (j=int(solV.size())-1 ; j >=0; j--)
		if (solV[j].maxR > temp) solV.erase(solV.begin() + j);
		for (j = int(solV.size())-1; j >= 0; j--) {
		if (solV[j].maxR > temp) {
		solV.erase(solV.begin() + j);
		}
		}

		// then by average residual
		temp = solV[0].averR; // find the lowest average residual
		for (j = 1; j < int(solV.size()); j++) {
		if (solV[j].averR < temp) temp = solV[j].averR;
		if (solV[j].averR < temp) {
		temp = solV[j].averR;
		}
		}
		//toss out any solutins that have average resiudals highter than ceil(lowest average residual)
		// toss out any solutions that have average resiudals higher than ceil(lowest average residual)
		temp = ceil(temp);
		for (j=int(solV.size())-1 ; j >=0; j--)
		if (solV[j].averR > temp) solV.erase(solV.begin() + j);
		for (j = int(solV.size()) - 1; j >= 0; j--) {
		if (solV[j].averR > temp) {
		solV.erase(solV.begin() + j);
		}
		}

		// finally choose the transformation with the lowest theta angle as the final filter
		k = 0;
		for (j = 1; j < int(solV.size()); j++) {
		if (fabs(solV[j].theta) < fabs(solV[k].theta)) k=j;
		if (fabs(solV[j].theta) < fabs(solV[k].theta)) {
		k = j;
		}
		}

		// save solV[k]
		@@ -559,14 +586,18 @@ void IsisMain() {
		trans[i].limit[1] = solV[k].limit[1] / double(conFid.size());
		}// end of scan loop

		// translation from 8 to 8
		trans[7].theta = 0.0;
		trans[7].dx = 0.0;
		trans[7].dy = 0.0;

		//we now have seven transformations that convert from scan i to scan i+1 for i = 0 to 7
		//as scan 8 is the left most scan and has a nominally identity transfrom from its sample-line
		//sytem to the stiched cube sample-line system Comibine transformations so that each transforms
		// We now have seven transformations that convert from scan i to scan i+1 for i = 0 to 7
		// as scan 8 is the left most scan and has a nominally identity transform from its sample-line
		// system to the stitched cube sample-line system. Combine transformations so that each transforms
		// is from scan i to scan 8
		for (i = 0; i < 7; i++) {
		for (j = i + 1; j < 8; j++) {
		//make i equal to the combined transformation of i and j
		// Make i equal to the combined transformation of i and j
		trans[i].theta = atan2( cos(trans[j].theta)*sin(trans[i].theta) +
		sin(trans[j].theta)*cos(trans[i].theta),
		cos(trans[j].theta)*cos(trans[i].theta) -
		@@ -581,12 +612,9 @@ void IsisMain() {
		}
		}

		//translation from 8 to 8
		trans[7].theta = 0.0;
		trans[7].dx = 0.0;
		trans[7].dy = 0.0;

		//now lets find the extents of the stiched image

		// Now lets find the extents of the stitched image
		double minS = 1,
		maxS = panC[7]->sampleCount(),
		minL = 1,
		@@ -596,7 +624,7 @@ void IsisMain() {
		scanS = panC[i]->sampleCount();
		scanL = panC[i]->lineCount();

		//convert the four corner points to the scan 8 domain and determine the greatest extents
		// Convert the four corner points to the scan 8 domain and determine the greatest extents
		temp = cos(trans[i].theta) - sin(trans[i].theta) + trans[i].dx;
		if (temp < minS) minS = temp;
		if (temp > maxS) maxS = temp;
		@@ -624,23 +652,25 @@ void IsisMain() {
		if (temp > maxL) maxL = temp;
		}

		//update the transformations to make the minimum line = 1
		for (i=0; i<8; i++)
		// Update the transformations to make the minimum line = 1
		for (i = 0; i < 8; i++) {
		trans[i].dy += 1- minL;
		}

		maxL += ceil(1 - minL);
		minL = 1;

		//update the transformations to make the minimum sample = 1
		for (i=0; i<8; i++)
		// Update the transformations to make the minimum sample = 1
		for (i = 0; i < 8; i++) {
		trans[i].dx += 1- minS;
		}

		maxS += ceil(1 - minS);
		minS = 1;

		// trans[i].limit was previous calculated as the center of mass of the conjugate fiducials in the
		// scani system
		//transfrom it now into the stiched cube system to be used as a limit between scans
		// Transform it now into the stitched cube system to be used as a limit between scans
		for (i = 0; i < 7; i++) {
		temp = trans[i].limit[0]*cos(trans[i].theta) -
		trans[i].limit[1]*sin(trans[i].theta) + trans[i].dx;
		@@ -651,16 +681,18 @@ void IsisMain() {



		//finally shift and invert the transformations so that they covert from the rubber sheeted
		// sub-cubes back to the sub-scans
		//the result of this will be transformations that will convert from a set of coordinate systems
		// that differ only by an integral transformation back to sub-scans, each sub-scan will be
		// rubber sheeted into one of these systems and
		//the limit[0] values record where they belong in the stiched image/system
		// Finally shift and invert the transformations so that they convert from the rubber sheeted
		// sub-cubes back to the sub-scans. The result of this will be transformations that will convert
		// from a set of coordinate systems that differ only by an integral transformation back to
		// sub-scans, each sub-scan will be rubber-sheeted into one of these systems, and the limit[0]
		// values record where they belong in the stiched image/system.
		for (i = 0; i < 8; i++) {
		if (i==7)
		if (i == 7) {
		sampleFrom = 1.0;
		else sampleFrom = trans[i].limit[0];
		}
		else {
		sampleFrom = trans[i].limit[0];
		}
		//a shift so that the transform puts sampleFrom at sample 1 of the rubber sheeted sub cubes
		trans[i].dx += -sampleFrom+1;
		trans[i].theta = -trans[i].theta;
		@@ -669,11 +701,11 @@ void IsisMain() {
		trans[i].dx = temp;
		}

		//make the final maxes integral
		// Make the final maxes integral
		maxS = ceil(maxS);
		maxL = ceil(maxL);

		//attributes of input and ouput cubes to process class
		// Attributes of input and ouput cubes to process class
		CubeAttributeOutput att;
		CubeAttributeInput attI;

		@@ -683,10 +715,12 @@ void IsisMain() {
		att.setFileFormat( panC[0]->format() );
		att.setByteOrder( panC[0]->byteOrder() );
		att.setPixelType( panC[0]->pixelType() );
		if (panC[0]->labelsAttached())
		if (panC[0]->labelsAttached()) {
		att.setLabelAttachment(AttachedLabel);
		else
		}
		else {
		att.setLabelAttachment(DetachedLabel);
		}

		//define an output cube
		outputC.setDimensions(int(maxS), int(maxL), 1);
		@@ -707,15 +741,23 @@ void IsisMain() {
		scanL = panC[i]->lineCount();

		//define the sample range
		if (i==0) sampleTo = maxS;
		else sampleTo = trans[i-1].limit[0];
		if (i==7) sampleFrom = 1.0;
		else sampleFrom = trans[i].limit[0];
		if (i == 0) {
		sampleTo = maxS;
		}
		else {
		sampleTo = trans[i-1].limit[0];
		}
		if (i == 7) {
		sampleFrom = 1.0;
		}
		else {
		sampleFrom = trans[i].limit[0];
		}

		Interpolator bilinearInt(Interpolator::BiLinearType);
		ProcessRubberSheet rubberS;

		//use ProcessRubber sheet to create the sub-cube for this scan
		//use ProcessRubber sheet to create the sub-cube for this scanz
		Trans2d3p transform(trans[i].theta, trans[i].dx, trans[i].dy,
		int(sampleTo - sampleFrom),
		int(maxL));
		@@ -723,7 +765,8 @@ void IsisMain() {
		rubberS.SetOutputCube(tempFile.expanded(),
		att,
		transform.OutputSamples(),
		transform.OutputLines(),1);
		transform.OutputLines(),
		1);
		rubberS.Progress()->SetText(
		QObject::tr("Transforming Scan%1: ").arg(i+1).toStdString().c_str());
		rubberS.StartProcess(transform, bilinearInt);

isis/src/apollo/apps/apollopanstitcher/apollopanstitcher.xml

+8 −5

Original line number	Diff line number	Diff line
		@@ -33,6 +33,9 @@
		<change name="Lynn Weller" date="2012-01-22">
		Application category name changed from Import and Export to Apollo. Fixes mantis ticket #951.
		</change>
		<change name="Christopher Combs" date="2017-06-29">
		Moved initialization of trans[7] to before the loop that uses it. Fixes #4948.
		</change>
		</history>

		<category>

isis/src/base/objs/Process/Process.cpp

+1 −2

Original line number	Diff line number	Diff line
		@@ -775,4 +775,3 @@ namespace Isis {
		}

		} // end namespace isis