I need to calculate the angle subtended between the start and end point of a spiral and the entry tangent of the spiral. The spiral starts as a straight line and ends with a known radius.

As an example, you can use a spiral with a length of 156.025, and a spiral parameter of 395 (the radius can be calculated: Radius = Parameter^2/Length).

I can calculate an almost correct answer with the following code:

VB Code:
  1. Dim dTheta As Double           'Exit tangent of clothoide

  2. Dim dAlpha As Double           'Starting tangent

  3. Dim dX As Double               'Opposite side of triangle formed by Start/End points and entry tangent

  4. Dim dY As Double               'Adjacent side of triangle formed by Start/End points and entry tangent

  5. Dim dRadius As Double          'Clothoide radius

  6.  
  7. 'Calculate end radius of clothoide

  8. dRadius = dParameter ^ 2 / dLength

  9.  
  10. 'Calculate end tangent angle

  11. dTheta = dLength / (2 * dRadius)

  12.  
  13. 'Calculate sides of triangle

  14. dX = dLength - (dLength ^ 5 / (10 * ((2 * dLength * dRadius) ^ 2))) _

  15.              + (dLength ^ 9 / (216 * ((2 * dLength * dRadius) ^ 4))) _

  16.              - (dLength ^ 13 / (9360 * ((2 * dLength * dRadius) ^ 6))) _

  17.              + (dLength ^ 17 / (685440 * ((2 * dLength * dRadius) ^ 8)))

  18.  
  19. dY = (dLength ^ 3 / (3 * (2 * dLength * dRadius))) _

  20.      - (dLength ^ 7 / (42 * ((2 * dLength * dRadius) ^ 3))) _

  21.      + (dLength ^ 11 / (1320 * ((2 * dLength * dRadius) ^ 5))) _

  22.      - (dLength ^ 15 / (75600 * ((2 * dLength * dRadius) ^ 7))) _

  23.      + (dLength ^ 19 / (6894720 * ((2 * dLength * dRadius) ^ 9)))

  24.  
  25. 'Calculate missing angle

  26. dAlpha = Atn(dY / dX)

But this solution gives an answer with an accuracy of +/- 0.001 gon and I need one order of magnitude more (ie. +/- 0.0001 gon or 0.0000015 radians).

I have done a search on the forum and tried changing the doubles to variants and using the cdec() function according to various posts but I always get the same answer.

Can anybody see a mistake in the code or advise as to the best way to acheive the required accuracy in the calculation?

Many Thanks