
(*********************************************************************

	Copyright 1991 by Roman E. Maeder
	
	Adapted from
	Roman E. Maeder: Programming in Mathematica,
	Second Edition, Addison-Wesley, 1991.

	Permission is hereby granted to make copies of this file for
	any purpose other than direct profit, or as part of a
	commercial product, provided this copyright notice is left
	intact.  Sale, other than for the cost of media, is prohibited.
	
	Permission is hereby granted to reproduce part or all of
	this file, provided that the source is acknowledged.

 *********************************************************************)

	(* Enclosed with Mathematica by permission *)


(*:Version: Mathematica 2.0 *)

(*:Name: Graphics`ComplexMap.m` *)

(*:Title: CartesianMap and PolarMap *)

(*:Author:
	Roman E. Maeder
*)

(*:Keywords:
	Cartesian coordinates, polar coordinates, lines, mapping
*)

(*:Requirements: none. *)

(*:Warnings: none. *)

(*:Sources:
	Roman E. Maeder: Programming in Mathematica,
	Second Edition, Addison-Wesley, 1991.
*)

(*:Summary:
This package plots the images of Cartesian coordinate lines and
polar coordinate lines under a user-supplied function f.
*)

BeginPackage["Graphics`ComplexMap`"]

CartesianMap::usage = "CartesianMap[f, {x0, x1, (dx)}, {y0, y1,
(dy)}, options...] plots the image of the Cartesian coordinate
lines under the function f.  The default values of dx and dy are
chosen so that the number of lines is equal to the value of the
option PlotPoints of Plot3D[]."

PolarMap::usage = "PolarMap[f, {r0:0, r1, (dr)}, {phi0, phi1,
(dphi)}, options...] plots the image of the polar coordinate lines
under the function f.  The default values of dr and dphi are chosen
so that the number of lines is equal to the value of the option
PlotPoints of Plot3D[].  The default for the phi range is {0, 2Pi}."

Begin["`Private`"]

huge = 10.0^6

SplitLine[vl_] :=
	Module[{vll, pos, linelist = {}, low, high},
		vll = If[NumberQ[#], #, Indeterminate]& /@ vl;
		pos = Flatten[ Position[vll, Indeterminate] ];
		pos = Union[ pos, {0, Length[vll]+1} ];
		Do[ low = pos[[i]]+1;
		    high = pos[[i+1]]-1;
		    If[ low < high, AppendTo[linelist, Take[vll, {low, high}]] ],
		   {i, 1, Length[pos]-1}];
		linelist
	]

MakeLines[points_] :=
	Module[{lines, newpoints},
		newpoints = points /.
			{ z_?NumberQ :> huge z/Abs[z] /; Abs[z] > huge,
			  z_?NumberQ :> 0.0 /; Abs[z] < 1/huge,
			  DirectedInfinity[z_] :> huge z/Abs[z] };
		lines = Join[ newpoints, Transpose[newpoints] ];
		lines = Flatten[ SplitLine /@ lines, 1 ];
		lines = Map[ {Re[#], Im[#]}&, lines, {2} ];
		lines = Map[ Line, lines ];
		Graphics[lines]
	]

FilterOptions[ command_Symbol, opts___ ] :=
	Module[{keywords = First /@ Options[command]},
		Sequence @@ Select[ {opts}, MemberQ[keywords, First[#]]& ]
	]

CartesianMap[ f_, {x0_, x1_, dx_:Automatic}, {y0_, y1_, dy_:Automatic}, opts___ ] :=
	Module[ {x, y, points, plotpoints, ndx=N[dx], ndy=N[dy]},
		plotpoints = PlotPoints /. {opts} /. Options[Plot3D];
		If[ dx === Automatic, ndx = N[(x1-x0)/(plotpoints-1)] ];
		If[ dy === Automatic, ndy = N[(y1-y0)/(plotpoints-1)] ];
		points = Table[ N[f[x + I y]],
			{x, x0, x1, ndx}, {y, y0, y1, ndy} ];
		Show[ MakeLines[points], FilterOptions[Graphics, opts],
			AspectRatio->Automatic, Axes->Automatic ]
	] /; NumberQ[N[x0]] && NumberQ[N[x1]] && NumberQ[N[y0]] && NumberQ[N[y1]]
		 (NumberQ[N[dx]] || dx === Automatic) &&
		 (NumberQ[N[dy]] || dy === Automatic)

PolarMap[ f_, {r0_:0, r1_, dr_:Automatic}, {phi0_, phi1_, dphi_:Automatic}, opts___ ] :=
	Module[ {r, phi, points, plotpoints, ndr=dr, ndphi=dphi},
		plotpoints = PlotPoints /. {opts} /. Options[Plot3D];
		If[ dr === Automatic, ndr = N[(r1-r0)/(plotpoints-1)] ];
		If[ dphi === Automatic, ndphi = N[(phi1-phi0)/(plotpoints-1)] ];
		points = Table[ N[f[r Exp[I phi]]],
			{r, r0, r1, ndr}, {phi, phi0, phi1, ndphi} ];
		Show[ MakeLines[points], FilterOptions[Graphics, opts],
			AspectRatio->Automatic, Axes->Automatic ]
	] /; NumberQ[N[r0]] && NumberQ[N[r1]] && NumberQ[N[phi0]] && NumberQ[N[phi1]]
		 (NumberQ[N[dr]] || dr === Automatic) &&
		 (NumberQ[N[dphi]] || dphi === Automatic)

PolarMap[ f_, rRange_List, opts___Rule ] :=
	PolarMap[ f, rRange, {0, 2Pi}, opts ]

End[ ]   (* Graphics`ComplexMap`Private` *)

Protect[CartesianMap, PolarMap]

EndPackage[ ]   (* Graphics`ComplexMap` *)

(*:Limitations: none known. *)

(*:Tests:

*)

(*:Examples:

CartesianMap[ Cos,{0.2, Pi-0.2 },{-2,2}]

CartesianMap[ Cos, {0.2, Pi-0.2, (Pi-0.4)/19}, {-2,2,4/16}]

CartesianMap[ Exp, {-1,1,0.2},{-2,2,0.2}]

PolarMap[ Log, {0.1,10,0.5},{-3,3,0.15}]

PolarMap[ 1/Conjugate[#]&, {0.1,5.1,0.5},{-Pi,Pi,Pi/12}]

CartesianMap[ Zeta,{0.1,0.9},{0,20}]

SetOptions[ Plot3D, PlotPoints->25]; CartesianMap[ Zeta,{0.1,0.9},{0,20}]

PolarMap[ Sqrt,{1},{-Pi-0.0001,Pi}]

PolarMap[ Sqrt,{1},{-Pi+0.0001,Pi}]

PolarMap[ Identity,{1,2},{-Pi,Pi,Pi/12}]

PolarMap[ Identity,{1,2},{-Pi,Pi,Pi/12},
          AspectRatio -> 0.5, Axes -> None, Framed -> True ]

CartesianMap[ 1/Conjugate[#]&, {-2,2,4/19},{-2,2,4/19},PlotRange->All]

f[z_] := Indeterminate /; Abs[z] < 10.^-3
f[z_] := 10^100 /; Abs[z-2] < 10.^-3
f[z_] := DirectedInfinity[-2-I] /; Abs[z-(-2+I)] < 10.^-3
f[z_] := DirectedInfinity[] /; Abs[z-(2-2I)] < 10.^-3
f[z_] := z
CartesianMap[f,{-2,2,1/3},{-2,2,1/3},Axes->None]
CartesianMap[1/Conjugate[#]&,{-2,2,0.2},{-2,2,0.2}]

*)
