CloudCompare forum

Code: Select all

[09:31:42] [doActionFitQuadric] Quadric equation: z = 0.191049 + 3.54679e-05 * x + 0.00128381 * y + -0.00411162 * x^2 + 0.00332104 * x.y + -0.000473885 * y^2

Code: Select all

[09:31:28] [doActionFitQuadric] Quadric equation: z = 0.467229 + 3.54679e-05 * x + 0.00128381 * y + -0.00411162 * x^2 + 0.00332104 * x.y + -0.000473885 * y^2

Code: Select all

Set fso = CreateObject("Scripting.FileSystemObject")

Set objShell = CreateObject("Wscript.Shell")

strPath = objShell.CurrentDirectory &"\"

Set folder = fso.GetFolder(strPath)

Set files = folder.Files

COGX = 0.8245
COGY = 28.636
COGZ = 73.2375

For each folderIdx In files
    
	filename = folderIdx.Name
	name = fso.getbasename(filename)
	extension = fso.GetExtensionName(filename)
	reconstructed = name & "_reconstructed.asc"

	If ( StrComp(extension,"txt",vbTextCompare) = 0 and not fso.FileExists(reconstructed) ) Then

		Set listFile = fso.OpenTextFile(filename)
			
		Set outputFile = fso.CreateTextFile(reconstructed, True)
		
		
			do while not listFile.AtEndOfStream 

				coordinates = split(listFile.ReadLine())
				
				newZ = Quadric(coordinates(0) - COGX, coordinates(1) - COGY) + COGZ
				
				outputFile.WriteLine(coordinates(0) & " " & coordinates(1) & " " & newZ)
			loop
		
			outputFile.Close
		
		End If
	
Next

Function Quadric(ByVal x, ByVal y)
	Dim z

	z = 0.467229 + 3.54679e-05 * x + 0.00128381 * y + -0.00411162 * x^2 + 0.00332104 * x * y + -0.000473885 * y^2
 
 	Quadric = z
End Function

Code: Select all

//the 4x4 transformation displayed in the Console
const ccGLMatrix& trans = quadric->getTransformation();
//the inverse transformation
ccGLMatrix invTrans = trans.inverse();

//for each point
for (unsigned i=0; i<newCloud->size(); ++i)
{
	//original point
	CCVector3* P = const_cast<CCVector3*>(newCloud->getPoint(i));
	//point in the local coordinate system
	CCVector3 Q = invTrans * (*P);
	//project the point on the quadric
	Q.z = eq[0] + eq[1]*Q.x + eq[2]*Q.y + eq[3]*Q.x*Q.x + eq[4]*Q.x*Q.y + eq[5]*Q.y*Q.y;
	//back project to the original coordinate system
	*P = trans * Q;
}