%------------------------------------------------------%  
%---Helix curve (parametrized on arc lenght)-----------% 
%---and intrinsic vectors------------------------------%  
%------------------------------------------------------%  
DEF helix (a,b::IsReal) = <K:a * cos ~ u, K:a * sin ~ u, K:b * u>
WHERE 
  u = ID / SQRT ~ K:(a*a + b*b)
END;

DEF helixtangent (a,b::IsReal) = AA:D:(helix:<a,b>);
DEF helixnormal (a,b::IsReal) = 
    K:((a*a + b*b)/a) scalarVectProd AA:(D ~ D):(helix:<a,b>);
DEF helixbinormal (a,b::IsReal) = helixtangent:<a,b> VectProd helixnormal:<a,b>;

DEF helixcurvatureCenter (a,b::IsReal) = 
  helix:<a,b> vectSum (K:((sqr:a + sqr:b)/a) scalarVectProd helixnormal:<a,b>) ;

DEF helixOsculCircle (a,b,s::IsReal) = 
  (T:<1,2,3>:center ~ Rotn:<angle, axis>):circle
WHERE
  circle = MAP:((circle3D ~ radius):<a,b, s>):(intervals:(2*PI):24),
  angle  = (- ~ ACOS ~ InnerProd):<ortho,<0,0,1>>,
  axis   = CONS:(helixnormal:<a,b>):s,
  center = CONS:(helixCurvatureCenter:<a,b>):s,
  ortho  = CONS:(helixbinormal:<a,b>):s
END;
DEF circle3D (r::IsReal) = [K:r * COS, K:r * SIN, K:0] ~ s1;
DEF radius (a,b,s::IsReal) = (VectNorm ~ VectDiff):
  < CONS:(helix:<a,b>):s, CONS:(helixcurvatureCenter:<a,b>):s >;
DEF CurveGraph (f::IsSeqOf:IsFun) = MAP:(CONS:f ~ s1):(intervals:(6*PI):90);

DEF assembly = STRUCT:<
  (CurveGraph ~ helix):<1,2>,
  (CurveGraph ~ helixcurvatureCenter):<1,2>,
  (JOIN ~ helixOsculCircle):<1,2, 2*PI>
>;
DEF out = (STRUCT ~ [ID, BOX:<1,2>]):assembly