BeginPackage["Caustic`WaveFront`"]

(*  File name: WaveFront.m                                 *)

(*  This package WaveFront.m was used for the              *)
(*  illustrations for my paper  O. Knill:                  *)
(*  "On nonconvex caustics of convex billiards"            *)
(*  which appeared in "Elemente der Mathematik, 1998       *)
(*  The package is not tuned for efficiency and clarity    *)
(*  It is not documented.                                  *)

(*  The first Version was developed at Caltech. The        *)
(*  second Version at the university of Arizona            *)

(*  Oliver Knill                                           *)
(*  Department of Mathematics                              *)
(*  University of Texas, Austin, TX                        *)
(*  e-mail:  knill@math.utexas.edu                         *)
(*  http://www.ma.utexas.edu/~knill                        *)

Figure12a              ::usage="";
Figure12b              ::usage="";
Figure14               ::usage="";

Begin["`Private`"] 

<<"Graphics`SurfaceOfRevolution`"

p=4; c={0,0,1}; P=N[Pi]; 
TR[t_]:=-I*p*(Exp[I t]-1) -
   I*Sum[          c[[k]] *(1/( k+1))*(Exp[I( k+1)*t]-1),{k,2,Length[c]}] -
   I*Sum[Conjugate[c[[k]]]*(1/(-k+1))*(Exp[I(-k+1)*t]-1),{k,2,Length[c]}];
f[alpha_]:=Re[TR[alpha]]; g[alpha_]:=Im[TR[alpha]];
T[alpha_,phi_]:={Cos[phi] f[alpha], -Sin[phi] f[alpha], g[alpha]};   
n[gamma_,delta_]:=-Sign[p]*(T[gamma+P,delta]-T[gamma,delta])/6;

PP:=SurfaceOfRevolution[{f[u],g[u]},{u,0,2 Pi},
  PlotPoints->20,AspectRatio->1]; 
SR[k_]:=Show[PP,ViewPoint->{5,5,k},
   Boxed->False,Axes->False,PlotRange->{{-p,p},{-p,p},{0,2.5p}}]; 

Tisch[t_]:={Re[TR[t]],Im[TR[t]]};
Tisch3D[t_]:={Re[TR[t]],0,Im[TR[t]]};
Tisch3D[t_,a_]:={{ Cos[a],Sin[a],0},
                 {-Sin[a],Cos[a],0},
                      {0,      0,1}}.{Re[TR[t]],0,Im[TR[t]]};
ShowIt:=Show[%,DisplayFunction->$DisplayFunction]; 

ThreeDTisch[m_]:=ParametricPlot3D[T[t,a],{t,Pi,2Pi},{a,0,2Pi},
  PlotPoints->m,Boxed->False,Axes->False,DisplayFunction->Identity];
ThreeDTisch1[m_]:=ParametricPlot3D[Tisch3D[t,a],{t,Pi,2Pi},{a,0,2Pi},
  PlotPoints->m,Boxed->False,Axes->False,DisplayFunction->Identity,
  FaceGrids->None];

ShowNormals[nn_]:=Module[{},
  s=Table[Line[{T[P k/nn,2P l/nn],T[P k/nn,2P l/nn]+n[P k/nn,2P l/nn]}],
  {k,nn},{l,nn}];
  Show[Graphics3D[Flatten[s]], 
    Boxed->False, 
    Axes->False, 
    AspectRatio->1,
    DisplayFunction->Identity]];

p=3; 
w[t_]:= N[I*Sum[(k/(k+1))*c[[k]]            (Exp[I( k+1)t]-1),{k,2,Length[c]}]+
        I*Sum[(-k/(-k+1))*Conjugate[c[[k]]] (Exp[I(-k+1)t]-1),{k,2,Length[c]}]+
        I*(p+2*Sum[c[[k]],{k,2,Length[c]}])];
Caustic3D[cc_List,m_]:=Module[{},c=cc;
  CC[t_,a_]:={{ Cos[a],Sin[a],0},
              {-Sin[a],Cos[a],0},
                    {0,     0,1}}.{Re[w[t]],0,Im[w[t]]};
  ParametricPlot3D[CC[t,a],{t,Pi,2Pi},{a,0,2Pi},
  PlotPoints->m,
  Boxed->False,Axes->False,
  AspectRatio->1]];

(* Plot the dimensional wave fronts  *)

Front[pp_]:=Module[{},p=pp; (* Print["newfront"]; *) ThreeDTisch[30]];
Front1[pp_]:=Module[{},p=pp;       a1=ThreeDTisch[10]; 
                       p=pp+0.001; a2=ThreeDTisch1[10]; {a1,a2}]; 
Front2[pp_]:=Module[{},p=pp; Print["newfront"]; 
  Show[{ShowNormals[20],ThreeDTisch[20]}]];
                              
Figure12a:=Show[ThreeDTisch[30],
  DisplayFunction->$DisplayFunction,
  ViewPoint->{5,5,-3}];
Figure12b:=Show[Caustic3D[{0,0,1},20], 
  DisplayFunction->$DisplayFunction,
  ViewPoint->{5,5,-3}]; 
Figure14:=Show[GraphicsArray[Partition[Table[Front2[1.5 k/8],{k,-8,7}],4]],
  DisplayFunction->$DisplayFunction]; 

End[]             
EndPackage[]

