(*--------------------------------------------------------------------------*)
(*----------------------------- FRACTALS -----------------------------------*)
(*--------------------------------------------------------------------------*)

(* All fractal curve approximations are treated by this package as lists of *)
(* segments, a segment being a list of two points, each of them being       *)
(* a list of two coordinates in the plane. The package contains a few       *)
(* examples of simple lists of segments (patterns). By replacing each of    *)
(* the segments of a list by a pattern, a new list is obtained, which can   *) 
(* be used as another pattern, or a new starting point. This package con-   *)
(* tains the replacement and iteration functions that permit the production *)
(* of many different curves, and also their representation.                 *)
(*--------------------------------------------------------------------------*)

BeginPackage["UVW`Fractals`"]

(*--------------------- Examples of basic patterns -------------------------*)

Triangle::usage="
Triangle returns an equilateral triangle as a list of three segments in 
the plane."

Star::usage="
Returns a star as a list of six segments in the plane."

Island::usage="
Returns a notch and a triangle as a list of eight segments in the plane."

Battlement::usage="
Returns a battlement as a list of eight segments in the plane."

Hat::usage="
Hat[sharpness] returns a list of four segments in the plane. " 

Wy::usage="
Wy[angle1,length1,angle2,length2] returns a Y-shaped list of four segments 
in the plane. angle1,length1,angle2 and length2 are the parameters of the 
two branches." 

RSWy::usage="
RSWy[n] returns a random sample of size n of outputs of the function Wy."

RSTruncs::usage="
RSTruncs[n] returns a random sample of vertical segments in the plane."


(*--------- Functions for the replacement of segments by patterns ----------*)

TransformSegment::usage="
TransformSegment[segment,pattern] returns a list of segments obtained from 
pattern by moving it in the plane so as to make its ends fit with 
those of the initial segment."

IteratePattern::usage="
IteratePattern[listofsegments,pattern,n] applies TransformSegment to each 
of the segments in listofsegments. Iterates n times."

IterateRandomPatterns::usage="
IterateRandomPattern[listofsegments,patterns,probas,n] iterates
n times the following operation : for each of the segments in listofsegments
one of the patterns in patterns is chosen randomly according to a 
probability read in probas, and the segment is transformed accordingly
by TransformSegment." 

(*----------------------- Graphical representation -------------------------*)

DrawSegments::usage="
DrawSegments[listofsegments] plots the elements of listofsegments in the
plane."

(*--------------------------------------------------------------------------*)

Begin["`Private`"]

(*--------------------- Examples of basic patterns -------------------------*)

Triangle   = { { {1.,0.},{0.,0.} },
               { {0.,0.},{0.5,N[Sqrt[3]/2]} }, 
               { {0.5,N[Sqrt[3]/2]},{1.,0.} } };

Island     = { { {0.,0.},{1.,0.} },
               { {1.,0.},{2.,-N[Sqrt[3]/2]} }, 
               { {2.,-N[Sqrt[3]/2]},{3.,-N[Sqrt[3]/2]} }, 
               { {2.,0.},{2.5,N[Sqrt[3]/2]} },
               { {2.5,N[Sqrt[3]/2]},{3.,0.} },
               { {3.,0.},{2.,0.} },
               { {3.,-N[Sqrt[3]/2]},{4.,0.} },
               { {4.,0.},{5.,0.} } };

Star       = { { {0.,0.},{1.,0.} },
               { {0.,0.},{0.5,N[Sqrt[3]/2]} },  
               { {0.,0.},{-0.5,N[Sqrt[3]/2]} },
               { {0.,0.},{-1.,0.} },
               { {0.,0.},{-0.5,-N[Sqrt[3]/2]} },  
               { {0.,0.},{0.5,-N[Sqrt[3]/2]} } };

Battlement = { { {0.,0.},{1.,0.} },
               { {1.,0.},{1.,1.} },  
               { {1.,1.},{2.,1.} },
               { {2.,1.},{2.,0.} },
               { {2.,0.},{2.,-1.} },  
               { {2.,-1.},{3.,-1.} }, 
               { {3.,-1.},{3.,0.} },
               { {3.,0.},{4.,0.} } };  

Hat[sharpness_]:= 
             { { {0.,0.},{1.,0.} },
               { {1.,0.},{1.5,N[sharpness]} },
               { {1.5,N[sharpness]},{2.,0.} }, 
               { {2.,0.},{3.,0.} } };

Wy[angle1_,length1_,angle2_,length2_]:=
             { { {0.,0.},{0.,1.} },
               { {0.,1.},{N[length1*Sin[angle1]],
                           1.+N[length1*Cos[angle1]]} },
               { {0.,1.},{-N[length2*Sin[angle2]],
                           1.+N[length2*Cos[angle2]]} },
               { {0.,1.3},{0.,2.} } };

RSWy[n_Integer]:=Table[Wy[N[Pi/2]*Random[],Random[],
                          N[Pi/2]*Random[],Random[]],{n}];

RSTruncs[n_Integer] := 
    Block[{abs,ord,height,truncs},

          truncs={};
          Do[(
             abs     = Random[Real,10.]; 
             ord     = Random[Real,3.];
             height  = Random[Real,{4.,7.}];
             truncs  = Append[truncs,{{abs,ord},{abs,ord+height}}];
             ),{n}];

          Return[truncs];
         ];

(*--------- Functions for the replacement of segments by patterns ----------*)

TransformSegment[segment_List,pattern_List]:=

    Block[{i,firstofpat,lastofpat,normofpat,
             firstofseg,lastofseg,normofseg,
             sinangle,cosangle,
             point, transpat},

              (* First and last point met in pattern and their distance-----*)

          firstofpat = First[Flatten[pattern,1]];
          lastofpat  = Last[Flatten[pattern,1]];
          normofpat  = N[Sqrt[ (lastofpat-firstofpat).
                               (lastofpat-firstofpat) ]];
                              
              (* End points and length of segment---------------------------*)

          firstofseg = segment[[1]];
          lastofseg  = segment[[2]];
          normofseg  = N[Sqrt[ (lastofseg-firstofseg).
                               (lastofseg-firstofseg)]];

              (* Sine and Cosine of the angle between segment and pattern---*)

          cosangle   = ((lastofseg-firstofseg).(lastofpat-firstofpat))
                              /(normofseg*normofpat);
          sinangle   = (({{0.,-1.},{1.,0.}}.(lastofseg-firstofseg))
                       .(lastofpat-firstofpat))
                              /(normofseg*normofpat);
               
              (* Transformations that map the end points of pattern---------*)
              (* onto those of segment--------------------------------------*)

          translation[point_List]:=point-(firstofpat-firstofseg);

          rotation[point_List]:=({{cosangle,sinangle},
                                  {-sinangle,cosangle}}
                                 .(point-firstofseg)) 
                                 + firstofseg;
                                  
          homothety[point_List]:=(point-firstofseg)
                                   *normofseg/normofpat + firstofseg;

              (* Transformation of pattern point by point-------------------*)
          transpat=Flatten[pattern,1];
          transpat=Table[homothety[rotation[translation[transpat[[j]] ]]],
                         {j,1,Length[transpat]}];
          transpat=Partition[transpat,2];
         
          Return[transpat];
         ];   (* End of function TransformSegment---------------------------*)


IteratePattern[listofsegments_List,pattern_List,n_Integer]:=

    Block[{i,jagged},

          jagged = listofsegments;

          Do[(
              (* Replace all segments in jagged by a transformed pattern----*)

              jagged = Table[TransformSegment[jagged[[i]],pattern],
                             {i,1,Length[jagged]}];
              jagged = Flatten[jagged,1];
             ),{n}];

          Return[jagged];
         ];   (* End of function IteratePattern-----------------------------*)


IterateRandomPatterns[listofsegments_List,patterns_List,
                                          probas_List,n_Integer]:=

    Block[{i,j,cumul,jagged},

              (* Distribution function of probas----------------------------*)
          cumul = N[probas];
          cumul = cumul/Apply[Plus,cumul];
          cumul = Rest[FoldList[Plus,0,cumul]];
          
              (* Definition of random choice according to probas------------*)
          choice := 
            Block[{j,ran},

                   ran=Random[];
                   j=1;
                   While[ ran>cumul[[j]] , j=j+1];
                   Return[j];
                 ];
                 
          jagged = listofsegments;

          Do[(
              (* Replace all segments in jagged by a transformed pattern----*)

             jagged = Table[TransformSegment[jagged[[i]],
                                             patterns[[choice]]
                                            ],
                            {i,1,Length[jagged]}];
             jagged = Flatten[jagged,1];
             ),{n}];

          Return[jagged];
         ];   (* End of function IteratePattern-----------------------------*)


(*----------------------- Graphical representation -------------------------*)

DrawSegments[listofsegments_List,opts___]:=
    Block[{i},

          Show[Graphics[Table[Line[listofsegments[[i]]],
                        {i,1,Length[listofsegments]}]], opts,
               Frame->True, FrameTicks->None]
         ]    (* End of function DrawSegments-------------------------------*)


End[ ]        (* End of Private Context-------------------------------------*)


EndPackage[ ]
(*--------------------------------------------------------------------------*)
