@stevebaer thank you for reaching me.
I will try to clear my problem I occurs once more.
What we are doing:
I work as a bridge engineer and we use parametric design for bridge design. Before we start project, we ask road designer for some inputs. Mostly it is alignment (3D road axis). However, most of the civil softwares work with the LandXML (.xml) format which is easily transferable. Therefore, we developed our .gha library which is capable to import LandXML file and create a very accurate 3D curve (alignment). Bridge shape is derived from that alignment, therefore we need very accurate and smooth curve.
How LandXML file is imported:
Alignment in LandXML file is written as a two lists of mathematical definitions (circular arcs, transition curves, straight line etc.) for road profile (vertical profile of road axis) and road axis in plan (projection of road axis to XY plan). Therefore, file is read, each subcurve in both lists is “recreated” from mathematical definition. Result of this routine are two separated curves (vertical profile and plan axis).
First curve is planar (vertical profile) and located in XZ plane near origin.
Second curve (plan projection) is planar as well and is located in XY plane in world coordinates according to coordination system (relatively in long distance from origin).
Final step - creation of final 3D curve
In this stage we need to create a final single one curve (3D curve) from both curves. It is done by flowing vertical profile over plan projection. This step is done by following routine (we do not use FlowSpaceMorph Class ):
- Length of vertical profile curve’s planar projection is calculated
- Based on a step value (component’s input) points in the X axis is created. These points represents points in different stationing (in civil engineering stationing in road is calculated according to 2D plan projection not real 3D road axis).
- From these points vertical (in the direction if Z axis) lines are created. Each lines intersect profile curve at some point.
- Intersection of the profile curve and each vertical line is calculated (with CurveCurve or CurveLine intersection method). From obtained intersection event we retrieve a Z coordinate of intersection point. At this point we have a elevations of the road axis in defined stationings.
- We calculate points at defined stationing from second curve (plan projection) with and obtain X and Y coordinates. This is done with Curve.PointAt Method
- Together with Z coodinate from profile curve we have a list of points which defines final 3D curve. This curve is done with Curve.CreateInterpolatedCurve Method
Problem
This approach worked so far. But for the first time we observed that Z coodinate in fourth step is calculated incorrectly. After investigation I have found that CurveCurve or CurveLine “fails” to correctly calculate Z coordinate of intersection point. If this approach is done with native components in Grasshopper Z coordinates are correct.
Here is part of my C# code with flowing routine. I am not a pro-coder so code can be little unreadable in some parts.
private SortedList<int,Point3d> FlowCurve(Curve Original, Curve Final_shape)
{
if (CancellationToken.IsCancellationRequested) { return new SortedList<int, Point3d>(); }
//Transform X_project = Transform.PlanarProjection(new Plane(new Point3d(0, 0, 0), new Vector3d(0, 0, 1)));
//Curve Vertical_alignment_projection = Original.DuplicateCurve();
//Vertical_alignment_projection.Transform(X_project);
Curve Vertical_alignment_projection = new Line(new Point3d(Original.PointAtStart.X, Original.PointAtStart.Y, Original.PointAtStart.Y), new Point3d(Original.PointAtEnd.X, Original.PointAtEnd.Y, Original.PointAtEnd.Y)).ToNurbsCurve();
SortedList<int, Point3d> pntsDict = new SortedList<int, Point3d>();
double min = Math.Min(Vertical_alignment_projection.GetLength(), Final_shape.GetLength() - (sectionStartStationing - planStartStationing));
double r = min / step;
int range = Convert.ToInt32(Math.Floor(r))+1;
List<int> finalRange = Enumerable.Range(0, range).ToList();
(onwerComponet as XMLGHCreateAlignment).UpdateMessage("Flowing curves - " + "\n" + "this may take a while.");
// Ošetření, kdyby jedna z křivek končila dříve
if ((onwerComponet as XMLGHCreateAlignment).MultiThreading)
{
Parallel.ForEach(finalRange, R =>
{
// Tady to funguje zvláštně.
if (CancellationToken.IsCancellationRequested) { return; }
Curve VerAligProjection = Original.DuplicateCurve();
Curve Finalshape = Final_shape.DuplicateCurve();
Point3d Pnt1 = VerAligProjection.PointAtLength(R * step);
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Line Vertical_line = new Line(Pnt1, Pnt2);
Point3d Pnt_horizontal_alignment = default;
Point3d Pnt_vertical_alignment = Intersection.CurveLine(Original, Vertical_line, Rhino.RhinoMath.ZeroTolerance, Rhino.RhinoMath.ZeroTolerance)[0].PointA;
//Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
if (Math.Round((sectionStartStationing - planStartStationing), 6) >= 0 || sectionStartStationing - planStartStationing + R * step >= 0)
{
if ((R == 0) && Math.Round((sectionStartStationing - planStartStationing), 6) == 0.0)
{
Pnt_horizontal_alignment = Finalshape.PointAtStart;
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
}
else
{
Pnt_horizontal_alignment = Finalshape.PointAtLength((sectionStartStationing - planStartStationing) + R * step);
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
}
lock (pntsDict)
{
if (CancellationToken.IsCancellationRequested) { return; }
pntsDict.Add(R, Pnt_horizontal_alignment);
};
}
});
}
else
{
foreach (int R in finalRange)
{
if (CancellationToken.IsCancellationRequested) { return new SortedList<int, Point3d>(); }
Curve VerAligProjection = Original.DuplicateCurve();
Curve Finalshape = Final_shape.DuplicateCurve();
Point3d Pnt1 = VerAligProjection.PointAtLength(R * step);
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Point3d Pnt_horizontal_alignment = default;
Line Vertical_line = new Line(Pnt1, Pnt2);
// Tohle vypočítá špatnou hodnotu. Zkusit nějak ještě zpřesnit pomocí Curve.ClosestPoint u vypočítaného a případěn křivku loftnout
// a udělat s Brep..... (nutno vyzkoušet a očůrat dané řešení)
// Do výstupu a popisu napsat problém a případně do výstupu dát půdorys a niveletu..
Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
if (Math.Round((sectionStartStationing - planStartStationing), 6) >= 0 || sectionStartStationing - planStartStationing + R * step >= 0)
{
if ((R == 0) && Math.Round((sectionStartStationing - planStartStationing), 6) == 0.0)
{
Pnt_horizontal_alignment = Finalshape.PointAtStart;
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
}
else
{
Pnt_horizontal_alignment = Finalshape.PointAtLength((sectionStartStationing - planStartStationing) + R * step);
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
}
lock (pntsDict)
{
if (CancellationToken.IsCancellationRequested) { return new SortedList<int, Point3d>(); }
pntsDict.Add(R, Pnt_horizontal_alignment);
};
}
}
}
if (CancellationToken.IsCancellationRequested) {return new SortedList<int, Point3d>(); }
double Len = default;
// Řešení posledního bodu
if ((sectionStartStationing - planStartStationing) >= 0)
{
// Podélňák začíná pozdějí a končí dříve
if (Vertical_alignment_projection.GetLength() + (sectionStartStationing - planStartStationing) < Final_shape.GetLength())
{
Len = Vertical_alignment_projection.GetLength();
Point3d Pnt1 = Vertical_alignment_projection.PointAtLength(Len);
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Line Vertical_line = new Line(Pnt1, Pnt2);
Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
Point3d Pnt_horizontal_alignment = Final_shape.PointAtLength(Len + (sectionStartStationing - planStartStationing));
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
pntsDict.Add(range, Pnt_horizontal_alignment);
}
// Podélňák začíná pozdějí a končí pozdějí
else
{
Len = Final_shape.GetLength();
Point3d Pnt1 = Vertical_alignment_projection.PointAtLength(Len - (sectionStartStationing - planStartStationing));
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Line Vertical_line = new Line(Pnt1, Pnt2);
Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
Point3d Pnt_horizontal_alignment = Final_shape.PointAtLength(Len);
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
pntsDict.Add(range, Pnt_horizontal_alignment); ;
}
}
else
{
// Podélňák začíná dříve a končí dříve
if (Vertical_alignment_projection.GetLength() + (sectionStartStationing - planStartStationing) < Final_shape.GetLength())
{
Len = Vertical_alignment_projection.GetLength();
Point3d Pnt1 = Vertical_alignment_projection.PointAtLength(Len);
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Line Vertical_line = new Line(Pnt1, Pnt2);
Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
Point3d Pnt_horizontal_alignment = Final_shape.PointAtLength(Len + (sectionStartStationing - planStartStationing));
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
pntsDict.Add(range, Pnt_horizontal_alignment);
}
// Podélňák začíná dříve a končí později
else
{
Len = Final_shape.GetLength();
Point3d Pnt1 = Vertical_alignment_projection.PointAtLength(Len - (sectionStartStationing - planStartStationing));
Point3d Pnt2 = new Point3d(Pnt1.X, Pnt1.Y, 10000);
Line Vertical_line = new Line(Pnt1, Pnt2);
Point3d Pnt_vertical_alignment = Intersection.CurveCurve(Original, Vertical_line.ToNurbsCurve(), AFRY_GH.Constants.Tolerances.DefaultTolerance, AFRY_GH.Constants.Tolerances.DefaultTolerance)[0].PointA;
Point3d Pnt_horizontal_alignment = Final_shape.PointAtLength(Len);
Pnt_horizontal_alignment.Z = Pnt_vertical_alignment.Z;
pntsDict.Add(range, Pnt_horizontal_alignment);
}
}
if (CancellationToken.IsCancellationRequested) { return new SortedList<int, Point3d>(); }
return pntsDict;
}
In the attachments you can find test file (.gha library with custom components, gh file with the importing workflow and problem description and xml file). Regarding the differences in Z coordinates, calculated with from the profile curve (native Grasshopper components) are generally correct (compared to the road design data and xml data). That means that difference is born when
CurveCurve or CurveLine method is used. However, in other cases this method worked fine.
I understand that difference cannot be precisely 0.0 but 0.011 is quite high. We usually work with mm so difference up to 0.001 is ok.
Sorry for the length of my post, but I feel important to describe my intentions.
Hope, my problem is clearer.
Thank you for your time.
Ondřej
Test_file.gh (33.9 KB)
GH_XML_Alignment.gha (50 KB)
SO101__Alignment.xml (12.7 KB)
