%D \module
%D   [       file=mp-tool.mpii,
%D        version=1998.02.15,
%D          title=\CONTEXT\ \METAPOST\ graphics,
%D       subtitle=auxiliary macros,
%D         author=Hans Hagen,
%D           date=\currentdate,
%D      copyright={PRAGMA ADE \& \CONTEXT\ Development Team}]
%C
%C This module is part of the \CONTEXT\ macro||package and is
%C therefore copyrighted by \PRAGMA. See mreadme.pdf for
%C details.

% def loadfile(expr name) = scantokens("input " & name & ";") enddef ;

if known context_tool : endinput ; fi ;

boolean context_tool ; context_tool := true ;

let @## = @# ;

%D New, version number testing:
%D
%D \starttyping
%D fill fullcircle scaled 2cm withcolor if mpversiongt("0.6") : red  else : green fi ;
%D fill fullcircle scaled 1cm withcolor if mpversionlt(0.6)   : blue else : white fi ;
%D \stoptyping

if not known mpversion : string mpversion ; mpversion := "0.641" ; fi ;

% newinternal metapostversion ; metapostversion := scantokens(mpversion) ;

newinternal metapostversion ; metapostversion := 2.0 ;

% vardef mpversiongt(expr s) =
%     scantokens (mpversion & " > " & if numeric s : decimal s else : s fi)
% enddef ;
% vardef mpversionlt(expr s) =
%     scantokens (mpversion & " < " & if numeric s : decimal s else : s fi)
% enddef ;
% vardef mpversioneq(expr s) =
%     scantokens (mpversion & " = " & if numeric s : decimal s else : s fi)
% enddef ;

%D More interesting:
%D
%D \starttyping
%D fill fullcircle scaled 4cm withcolor if mpversiongt("0.6")     : red    else : green fi ;
%D fill fullcircle scaled 2cm withcolor if mpversionlt(0.6)       : blue   else : white fi ;
%D fill fullcircle scaled 1cm withcolor if mpversioncmp(0.6,">=") : yellow else : black fi ;
%D \stoptyping

vardef mpversioncmp(expr s, c) =
    scantokens (mpversion & c & if numeric s : decimal s else : s fi)
enddef ;

vardef mpversionlt (expr s) = mpversioncmp(s, "<") enddef ;
vardef mpversioneq (expr s) = mpversioncmp(s, "=") enddef ;
vardef mpversiongt (expr s) = mpversioncmp(s, ">") enddef ;

%D We always want \EPS\ conforming output, so we say:

prologues    := 1 ;
warningcheck := 0 ;
mpprocset    := 1 ;

%D Namespace handling:

% let exclamationmark = ! ;
% let questionmark    = ? ;
%
% def unprotect =
%   let ! = relax ;
%   let ? = relax ;
% enddef ;
%
% def protect =
%   let ! = exclamationmark ;
%   let ? = questionmark ;
% enddef ;
%
% unprotect ;
%
% mp!some!module = 10 ; show mp!some!module ; show somemodule ;
%
% protect ;

string space ; space := char 32 ;
string CRLF  ; CRLF  := char 10 & char 13 ;

vardef ddecimal primary p =
    decimal xpart p & " " & decimal ypart p
enddef ;

%D Plain compatibility:

string plain_compatibility_data ; plain_compatibility_data := "" ;

def startplaincompatibility =
    begingroup ;
    scantokens plain_compatibility_data ;
enddef ;

def stopplaincompatibility =
    endgroup ;
enddef ;

% is now built in
%
% extra_endfig := extra_endfig
%   & "special "
%        & "("
%        & ditto
%        & "%%HiResBoundingBox: "
%        & ditto
%        & "&ddecimal llcorner currentpicture"
%        & "&space"
%        & "&ddecimal urcorner currentpicture"
%        & ");";

%D More neutral:

let triplet    = rgbcolor ;
let quadruplet = cmykcolor ;

%D Crap (experimental, not used):

def forcemultipass =
    % extra_endfig := extra_endfig & "special(" & ditto  & "%%MetaPostOption: multipass" & ditto & ");" ;
enddef ;

%D Colors:

let grayscale = numeric ;
let greyscale = numeric ;

vardef colorpart expr c =
    if not picture c :
        0
    elseif colormodel c = greycolormodel :
        greypart c
    elseif colormodel c = rgbcolormodel  :
        (redpart c,greenpart c,bluepart c)
    elseif colormodel c = cmykcolormodel :
        (cyanpart c,magentapart c,yellowpart c,blackpart c)
    else :
        0 % black
    fi
enddef ;

vardef colorlike(text c) text v = % colorlike(a) b, c, d ;
    save _p_ ; picture _p_ ;
    forsuffixes i=v :
        _p_ := image(draw origin withcolor c ;) ; % intercept pre and postscripts
        if (colormodel _p_ = cmykcolormodel) :
            cmykcolor i ;
        elseif (colormodel _p_ = rgbcolormodel) :
            rgbcolor i ;
        else :
            greycolor i ;
        fi ;
    endfor ;
enddef ;

%D Also handy (when we flush colors):

vardef dddecimal primary c =
    decimal redpart c & " " & decimal greenpart c & " " & decimal bluepart c
enddef ;

vardef ddddecimal primary c =
    decimal cyanpart c & " " & decimal magentapart c & " " & decimal yellowpart c & " " & decimal blackpart c
enddef ;

vardef colordecimals primary c =
    if cmykcolor c  :
        decimal cyanpart c & ":" & decimal magentapart c & ":" & decimal yellowpart c & ":" & decimal blackpart c
    elseif rgbcolor c :
        decimal redpart c & ":" & decimal greenpart c & ":" & decimal bluepart c
    else :
        decimal c
    fi
enddef ;

%D We have standardized data file names:

def job_name =
    jobname
enddef ;

def data_mpd_file =
    job_name & "-mp.mpd"
enddef ;

%D Because \METAPOST\ has a hard coded limit of 4~datafiles,
%D we need some trickery when we have multiple files.

if unknown collapse_data :
    boolean collapse_data ;
    collapse_data := false ;
fi ;

boolean savingdata     ; savingdata     := false ;
boolean savingdatadone ; savingdatadone := false ;

def savedata expr txt =
    write if collapse_data :
        txt
    else :
        if savingdata : txt else : "\MPdata{" & decimal charcode & "}{" & txt & "}" fi & "%"
    fi to data_mpd_file ;
enddef ;

def startsavingdata =
    savingdata := true ;
    savingdatadone := true ;
    if collapse_data :
        write "\MPdata{" & decimal charcode & "}{%" to data_mpd_file ;
    fi ;
enddef ;

def stopsavingdata =
    if collapse_data :
        write "}%" to data_mpd_file ;
    fi ;
    savingdata := false ;
enddef ;

def finishsavingdata =
    if savingdatadone :
        write EOF to data_mpd_file ;
        savingdatadone := false ;
    fi ;
enddef ;

%D Instead of a keystroke eating save and allocation
%D sequence, you can use the \citeer {new} alternatives to
%D save and allocate in one command.

def newcolor     text v = forsuffixes i=v : save i ; color     i ; endfor ; enddef ;
def newnumeric   text v = forsuffixes i=v : save i ; numeric   i ; endfor ; enddef ;
def newboolean   text v = forsuffixes i=v : save i ; boolean   i ; endfor ; enddef ;
def newtransform text v = forsuffixes i=v : save i ; transform i ; endfor ; enddef ;
def newpath      text v = forsuffixes i=v : save i ; path      i ; endfor ; enddef ;
def newpicture   text v = forsuffixes i=v : save i ; picture   i ; endfor ; enddef ;
def newstring    text v = forsuffixes i=v : save i ; string    i ; endfor ; enddef ;
def newpair      text v = forsuffixes i=v : save i ; pair      i ; endfor ; enddef ;

%D Sometimes we don't want parts of the graphics add to the
%D bounding box. One way of doing this is to save the bounding
%D box, draw the graphics that may not count, and restore the
%D bounding box.
%D
%D \starttyping
%D push_boundingbox currentpicture;
%D pop_boundingbox currentpicture;
%D \stoptyping
%D
%D The bounding box can be called with:
%D
%D \starttyping
%D boundingbox currentpicture
%D inner_boundingbox currentpicture
%D outer_boundingbox currentpicture
%D \stoptyping
%D
%D Especially the latter one can be of use when we include
%D the graphic in a document that is clipped to the bounding
%D box. In such occasions one can use:
%D
%D \starttyping
%D set_outer_boundingbox currentpicture;
%D \stoptyping
%D
%D Its counterpart is:
%D
%D \starttyping
%D set_inner_boundingbox p
%D \stoptyping

path    mfun_boundingbox_stack ;
numeric mfun_boundingbox_stack_depth ;

mfun_boundingbox_stack_depth := 0 ;

def pushboundingbox text p =
    mfun_boundingbox_stack_depth := mfun_boundingbox_stack_depth + 1 ;
    mfun_boundingbox_stack[mfun_boundingbox_stack_depth] := boundingbox p ;
enddef ;

def popboundingbox text p =
    setbounds p to mfun_boundingbox_stack[mfun_boundingbox_stack_depth] ;
    mfun_boundingbox_stack[mfun_boundingbox_stack_depth] := origin ;
    mfun_boundingbox_stack_depth := mfun_boundingbox_stack_depth - 1 ;
enddef ;

let push_boundingbox = pushboundingbox ; % downward compatible
let pop_boundingbox  = popboundingbox  ; % downward compatible

vardef boundingbox primary p =
    if (path p) or (picture p) :
        llcorner p -- lrcorner p -- urcorner p -- ulcorner  p
    else :
        origin
    fi -- cycle
enddef;

vardef innerboundingbox primary p =
    top  rt llcorner p --
    top lft lrcorner p --
    bot lft urcorner p --
    bot  rt ulcorner p -- cycle
enddef;

vardef outerboundingbox primary p =
    bot lft llcorner p --
    bot  rt lrcorner p --
    top  rt urcorner p --
    top lft ulcorner p -- cycle
enddef;

def inner_boundingbox = innerboundingbox enddef ;
def outer_boundingbox = outerboundingbox enddef ;

vardef set_inner_boundingbox text q = % obsolete
    setbounds q to innerboundingbox q;
enddef;

vardef set_outer_boundingbox text q = % obsolete
    setbounds q to outerboundingbox q;
enddef;

%D Some missing functions can be implemented rather straightforward (thanks to
%D Taco and others):

pi := 3.14159265358979323846 ; radian := 180/pi ; % 2pi*radian = 360 ;

% let +++ = ++ ;

numeric Pi ; Pi := pi ; % for some old compatibility reasons i guess

vardef sqr  primary x = x*x                                 enddef ;
vardef log  primary x = if x=0: 0 else: mlog(x)/mlog(10) fi enddef ;
vardef ln   primary x = if x=0: 0 else: mlog(x)/256 fi      enddef ;
vardef exp  primary x = (mexp 256)**x                       enddef ;
vardef inv  primary x = if x=0: 0 else: x**-1 fi            enddef ;

vardef pow (expr x,p) = x**p                                enddef ;

vardef tand   primary x = sind(x)/cosd(x)  enddef ;
vardef cotd   primary x = cosd(x)/sind(x)  enddef ;

vardef sin    primary x = sind(x*radian)   enddef ;
vardef cos    primary x = cosd(x*radian)   enddef ;
vardef tan    primary x = sin(x)/cos(x)    enddef ;
vardef cot    primary x = cos(x)/sin(x)    enddef ;

vardef asin   primary x = angle((1+-+x,x)) enddef ;
vardef acos   primary x = angle((x,1+-+x)) enddef ;
vardef atan   primary x = angle(1,x)       enddef ;

vardef invsin primary x = (asin(x))/radian enddef ;
vardef invcos primary x = (acos(x))/radian enddef ;
vardef invtan primary x = (atan(x))/radian enddef ;

vardef acosh  primary x = ln(x+(x+-+1))    enddef ;
vardef asinh  primary x = ln(x+(x++1))     enddef ;

vardef sinh   primary x = save xx ; xx = exp x ; (xx-1/xx)/2 enddef ;
vardef cosh   primary x = save xx ; xx = exp x ; (xx+1/xx)/2 enddef ;

%D Sometimes this is handy:

def undashed =
    dashed nullpicture
enddef ;

%D We provide two macros for drawing stripes across a shape.
%D The first method (with the n suffix) uses another method,
%D slower in calculation, but more efficient when drawn. The
%D first macro divides the sides into n equal parts. The
%D first argument specifies the way the lines are drawn, while
%D the second argument identifier the way the shape is to be
%D drawn.
%D
%D \starttyping
%D stripe_path_n
%D   (dashed evenly withcolor blue)
%D   (filldraw)
%D   fullcircle xscaled 100 yscaled 40 shifted (50,50) withpen pencircle scaled 4;
%D \stoptyping
%D
%D The a (or angle) alternative supports arbitrary angles and
%D is therefore more versatile.
%D
%D \starttyping
%D stripe_path_a
%D   (withpen pencircle scaled 2 withcolor red)
%D   (draw)
%D   fullcircle xscaled 100 yscaled 40 withcolor blue;
%D \stoptyping
%D
%D We have two alternatives, controlled by arguments or defaults (when arguments
%D are zero).
%D
%D The newer and nicer interface is used as follows (triggered by a question by Mari):
%D
%D \starttyping
%D draw image (draw fullcircle scaled 3cm shifted (0cm,0cm) withcolor green) numberstriped (1,10,3) withcolor red ;
%D draw image (draw fullcircle scaled 3cm shifted (3cm,0cm) withcolor green) numberstriped (2,20,3) withcolor green ;
%D draw image (draw fullcircle scaled 3cm shifted (3cm,3cm) withcolor green) numberstriped (3,10,5) withcolor blue ;
%D draw image (draw fullcircle scaled 3cm shifted (0cm,3cm) withcolor green) numberstriped (4,20,5) withcolor yellow ;
%D
%D draw image (draw fullcircle scaled 3cm shifted (6cm,0cm) withcolor green) anglestriped (1,20,2) withcolor red ;
%D draw image (draw fullcircle scaled 3cm shifted (9cm,0cm) withcolor green) anglestriped (2,40,2) withcolor green ;
%D draw image (draw fullcircle scaled 3cm shifted (9cm,3cm) withcolor green) anglestriped (3,60,2) withcolor blue ;
%D draw image (draw fullcircle scaled 3cm shifted (6cm,3cm) withcolor green) anglestriped (4,80,2) withcolor yellow ;
%D
%D draw image (
%D     draw fullcircle scaled 3cm shifted (0cm,0cm) withcolor green withpen pencircle scaled 2mm ;
%D     draw fullcircle scaled 2cm shifted (0cm,1cm) withcolor blue  withpen pencircle scaled 3mm ;
%D ) shifted (9cm,0cm) numberstriped (1,10,3) withcolor red ;
%D
%D draw image (
%D     draw fullcircle scaled 3cm shifted (0cm,0cm) withcolor green  withpen pencircle scaled 2mm ;
%D     draw fullcircle scaled 2cm shifted (0cm,1cm) withcolor blue   withpen pencircle scaled 3mm ;
%D ) shifted (12cm,0cm) numberstriped (2,10,3) withcolor red ;
%D
%D draw image (
%D     draw fullcircle scaled 3cm shifted (0cm,0cm) withcolor green  withpen pencircle scaled 2mm ;
%D     draw fullcircle scaled 2cm shifted (0cm,1cm) withcolor blue   withpen pencircle scaled 3mm ;
%D ) shifted (9cm,5cm) numberstriped (3,10,3) withcolor red ;
%D
%D draw image (
%D     draw fullcircle scaled 3cm shifted (0cm,0cm) withcolor green  withpen pencircle scaled 2mm ;
%D     draw fullcircle scaled 2cm shifted (0cm,1cm) withcolor blue   withpen pencircle scaled 3mm ;
%D ) shifted (12cm,5cm) numberstriped (4,10,3) withcolor red ;
%D \stoptyping

stripe_n     := 10;
stripe_slot  :=  3;
stripe_gap   :=  5;
stripe_angle := 45;

def mfun_tool_striped_number_action text extra =
    for i = 1/used_n step 1/used_n until 1 :
        draw point (1+i) of bounds -- point (3-i) of bounds withpen pencircle scaled penwidth extra ;
    endfor ;
    for i = 0 step 1/used_n until 1 :
        draw point (3+i) of bounds -- point (1-i) of bounds withpen pencircle scaled penwidth extra ;
    endfor ;
enddef ;

def mfun_tool_striped_set_options(expr option) =
    save isinner, swapped ;
    boolean isinner, swapped ;
    if option = 1 :
        isinner := false ;
        swapped := false ;
    elseif option = 2 :
        isinner := true ;
        swapped := false ;
    elseif option = 3 :
        isinner := false ;
        swapped := true ;
    elseif option = 4 :
        isinner := true ;
        swapped := true ;
    else :
        isinner := false ;
        swapped := false ;
    fi ;
enddef ;

vardef mfun_tool_striped_number(expr option, p, s_n, s_slot) text extra =
    image (
        begingroup ;
        save pattern, shape, bounds, penwidth, used_n, used_slot ;
        picture pattern, shape ; path bounds ; numeric used_s, used_slot ;
        mfun_tool_striped_set_options(option) ;
        used_slot := if s_slot = 0 : stripe_slot else : s_slot fi ;
        used_n := if s_n = 0 : stripe_n else : s_n fi ;
        shape := image(draw p) ;
        bounds := boundingbox shape ;
        penwidth := min(ypart urcorner shape - ypart llcorner shape, xpart urcorner shape - xpart llcorner shape) / (used_slot * used_n) ;
        pattern := image (
            if isinner :
                mfun_tool_striped_number_action extra ;
                for s within shape :
                    if stroked s or filled s :
                        clip currentpicture to pathpart s ;
                    fi
                endfor ;
            else :
                for s within shape :
                    if stroked s or filled s :
                        draw image (
                            mfun_tool_striped_number_action extra ;
                            clip currentpicture to pathpart s ;
                        ) ;
                    fi ;
                endfor ;
            fi ;
        ) ;
        if swapped :
            addto currentpicture also shape ;
            addto currentpicture also pattern ;
        else :
            addto currentpicture also pattern ;
            addto currentpicture also shape ;
        fi ;
        endgroup ;
    )
enddef ;

def mfun_tool_striped_angle_action text extra =
    for i = minimum -.5used_gap step used_gap until maximum :
        draw (minimum,i) -- (maximum,i) extra ;
    endfor ;
    currentpicture := currentpicture rotated used_angle ;
enddef ;

vardef mfun_tool_striped_angle(expr option, p, s_angle, s_gap) text extra =
    image (
        begingroup ;
        save pattern, shape, mask, maximum, minimum, centrum, used_angle, used_gap ;
        picture pattern, shape, mask ; numeric maximum, minimum ; pair centrum ; numeric used_angle, used_gap ;
        mfun_tool_striped_set_options(option) ;
        used_angle := if s_angle = 0 : stripe_angle else : s_angle fi ;
        used_gap := if s_gap = 0 : stripe_gap else : s_gap fi ;
        shape := image(draw p) ;
        centrum := center shape ;
        shape := shape shifted - centrum ;
        mask := shape rotated used_angle ;
        maximum := max (xpart llcorner mask, xpart urcorner mask, ypart llcorner mask, ypart urcorner mask) ;
        minimum := min (xpart llcorner mask, xpart urcorner mask, ypart llcorner mask, ypart urcorner mask) ;
        pattern := image (
            if isinner :
                mfun_tool_striped_angle_action extra ;
                for s within shape :
                    if stroked s or filled s :
                        clip currentpicture to pathpart s ;
                    fi
                endfor ;
            else :
                for s within shape :
                    if stroked s or filled s :
                        draw image (
                            mfun_tool_striped_angle_action extra ;
                            clip currentpicture to pathpart s ;
                        ) ;
                    fi ;
                endfor ;
            fi ;
        ) ;
        if swapped :
            addto currentpicture also shape ;
            addto currentpicture also pattern ;
        else :
            addto currentpicture also pattern ;
            addto currentpicture also shape ;
        fi ;
        currentpicture := currentpicture shifted centrum ;
        endgroup ;
    )
enddef;

newinternal striped_normal_inner  ; striped_normal_inner  := 1 ;
newinternal striped_reverse_inner ; striped_reverse_inner := 2 ;
newinternal striped_normal_outer  ; striped_normal_outer  := 3 ;
newinternal striped_reverse_outer ; striped_reverse_outer := 4 ;

secondarydef p anglestriped s =
    mfun_tool_striped_angle(redpart s,p,greenpart s,bluepart s)
enddef ;

secondarydef p numberstriped s =
    mfun_tool_striped_number(redpart s,p,greenpart s,bluepart s)
enddef ;

% for old times sake:

def stripe_path_n (text s_spec) (text s_draw) expr s_path =
    do_stripe_path_n (s_spec) (s_draw) (s_path)
enddef;

def do_stripe_path_n (text s_spec) (text s_draw) (expr s_path) text s_text =
    draw image(s_draw s_path s_text) numberstriped(3,0,0) s_spec ;
enddef ;

def stripe_path_a (text s_spec) (text s_draw) expr s_path =
    do_stripe_path_a (s_spec) (s_draw) (s_path)
enddef;

def do_stripe_path_a (text s_spec) (text s_draw) (expr s_path) text s_text =
    draw image(s_draw s_path s_text) anglestriped(3,0,0) s_spec ;
enddef ;

%D A few normalizing macros:
%D
%D \starttypen
%D xscale_currentpicture  ( width )
%D yscale_currentpicture  ( height )
%D xyscale_currentpicture ( width, height )
%D scale_currentpicture   ( width, height )
%D \stoptypen

% def xscale_currentpicture(expr the_width) =
%   natural_width  := xpart urcorner currentpicture - xpart llcorner currentpicture;
%   currentpicture := currentpicture scaled (the_width/natural_width) ;
% enddef;
%
% def yscale_currentpicture(expr the_height ) =
%   natural_height := ypart urcorner currentpicture - ypart llcorner currentpicture;
%   currentpicture := currentpicture scaled (the_height/natural_height) ;
% enddef;
%
% def xyscale_currentpicture(expr the_width, the_height) =
%   natural_width  := xpart urcorner currentpicture - xpart llcorner currentpicture;
%   natural_height := ypart urcorner currentpicture - ypart llcorner currentpicture;
%   currentpicture := currentpicture
%     xscaled (the_width/natural_width)
%     yscaled (the_height/natural_height) ;
% enddef;
%
% def scale_currentpicture(expr the_width, the_height) =
%   xscale_currentpicture(the_width) ;
%   yscale_currentpicture(the_height) ;
% enddef;

% nog eens uitbreiden zodat path en pic worden afgehandeld.

%   natural_width  := xpart urcorner currentpicture - xpart llcorner currentpicture;
%   currentpicture := currentpicture scaled (the_width/natural_width) ;

primarydef p xsized w =
    (p if (bbwidth (p)>0) and (w>0) : scaled (w/bbwidth (p)) fi)
enddef ;

primarydef p ysized h =
    (p if (bbheight(p)>0) and (h>0) : scaled (h/bbheight(p)) fi)
enddef ;

primarydef p xysized s =
    begingroup
    save wh, w, h ; pair wh ; numeric w, h ;
    wh := paired (s) ; w := bbwidth(p) ; h := bbheight(p) ;
    p
        if (w>0) and (h>0) :
            if xpart wh > 0 : xscaled (xpart wh/w) fi
            if ypart wh > 0 : yscaled (ypart wh/h) fi
        fi
    endgroup
enddef ;

let sized = xysized ;

def xscale_currentpicture(expr w) = % obsolete
    currentpicture := currentpicture xsized w ;
enddef;

def yscale_currentpicture(expr h) = % obsolete
    currentpicture := currentpicture ysized h ;
enddef;

def xyscale_currentpicture(expr w, h) = % obsolete
    currentpicture := currentpicture xysized (w,h) ;
enddef;

def scale_currentpicture(expr w, h) = % obsolete
    currentpicture := currentpicture xsized w ;
    currentpicture := currentpicture ysized h ;
enddef;

%D A full circle is centered at the origin, while a unitsquare
%D is located in the first quadrant. Now guess what kind of
%D path fullsquare and unitcircle do return.

path fullsquare, unitcircle ;

fullsquare := unitsquare shifted - center unitsquare ;
unitcircle := fullcircle shifted urcorner fullcircle ;

%D Some more paths:

path urcircle, ulcircle, llcircle, lrcircle ;

urcircle := origin -- (+.5,0) & (+.5,0){up}    .. (0,+.5) & (0,+.5) -- cycle ;
ulcircle := origin -- (0,+.5) & (0,+.5){left}  .. (-.5,0) & (-.5,0) -- cycle ;
llcircle := origin -- (-.5,0) & (-.5,0){down}  .. (0,-.5) & (0,-.5) -- cycle ;
lrcircle := origin -- (0,-.5) & (0,-.5){right} .. (+.5,0) & (+.5,0) -- cycle ;

path tcircle, bcircle, lcircle, rcircle ;

tcircle = origin -- (+.5,0) & (+.5,0) {up}    .. (0,+.5) .. {down}  (-.5,0) -- cycle ;
bcircle = origin -- (-.5,0) & (-.5,0) {down}  .. (0,-.5) .. {up}    (+.5,0) -- cycle ;
lcircle = origin -- (0,+.5) & (0,+.5) {left}  .. (-.5,0) .. {right} (0,-.5) -- cycle ;
rcircle = origin -- (0,-.5) & (0,-.5) {right} .. (+.5,0) .. {left}  (0,+.5) -- cycle ;

path urtriangle, ultriangle, lltriangle, lrtriangle ; % watch out: it's contrary to what you expect and starts in the origin

urtriangle := origin -- (+.5,0) -- (0,+.5) -- cycle ;
ultriangle := origin -- (0,+.5) -- (-.5,0) -- cycle ;
lltriangle := origin -- (-.5,0) -- (0,-.5) -- cycle ;
lrtriangle := origin -- (0,-.5) -- (+.5,0) -- cycle ;

path unitdiamond, fulldiamond ;

unitdiamond := (.5,0) -- (1,.5) -- (.5,1) -- (0,.5) -- cycle ;
fulldiamond := unitdiamond shifted - center unitdiamond ;

%D More robust:

% let  normalscaled =  scaled ;
% let normalxscaled = xscaled ;
% let normalyscaled = yscaled ;
%
% def  scaled expr s =  normalscaled (s) enddef ;
% def xscaled expr s = normalxscaled (s) enddef ;
% def yscaled expr s = normalyscaled (s) enddef ;

%D Shorter

primarydef p xyscaled q = % secundarydef does not work out well
    begingroup
    save qq ; pair qq ;
    qq = paired(q) ;
    p
%         if xpart qq <> 0 : xscaled (xpart qq) fi
%         if ypart qq <> 0 : yscaled (ypart qq) fi
        xscaled (xpart qq)
        yscaled (ypart qq)
    endgroup
enddef ;

%D Some personal code that might move to another module

def set_grid(expr w, h, nx, ny) =
    boolean grid[][] ; boolean grid_full ;
    numeric grid_w, grid_h, grid_nx, grid_ny, grid_x, grid_y, grid_left ;
    grid_w := w ;
    grid_h := h ;
    grid_nx := nx ;
    grid_ny := ny ;
    grid_x := round(w/grid_nx) ; % +.5) ;
    grid_y := round(h/grid_ny) ; % +.5) ;
    grid_left := (1+grid_x)*(1+grid_y) ;
    grid_full := false ;
    for i=0 upto grid_x :
        for j=0 upto grid_y :
            grid[i][j] := false ;
        endfor ;
    endfor ;
enddef ;

vardef new_on_grid(expr _dx_, _dy_) =
    dx := _dx_ ;
    dy := _dy_ ;
    ddx := min(round(dx/grid_nx),grid_x) ; % +.5),grid_x) ;
    ddy := min(round(dy/grid_ny),grid_y) ; % +.5),grid_y) ;
    if not grid_full and not grid[ddx][ddy] :
        grid[ddx][ddy] := true ;
        grid_left := grid_left-1 ;
        grid_full := (grid_left=0) ;
        true
    else :
        false
    fi
enddef ;

%D usage: \type{innerpath peepholed outerpath}.
%D
%D beginfig(1);
%D   def fullsquare = (unitsquare shifted -center unitsquare) enddef ;
%D   fill (fullsquare scaled 200) withcolor red ;
%D   path p ; p := (fullcircle scaled 100) ; bboxmargin := 0 ;
%D   fill p peepholed bbox p ;
%D endfig;

secondarydef p peepholed q =
    begingroup
    save start ; pair start ;
    start := point 0 of p ;
    if xpart start >= xpart center p :
        if ypart start >= ypart center p :
            urcorner q -- ulcorner q -- llcorner q -- lrcorner q --
            reverse  p -- lrcorner q -- cycle
        else :
            lrcorner q -- urcorner q -- ulcorner q -- llcorner q --
            reverse  p -- llcorner q -- cycle
        fi
    else :
        if ypart start > ypart center p :
            ulcorner q -- llcorner q -- lrcorner q -- urcorner q --
            reverse  p -- urcorner q -- cycle
        else :
            llcorner q -- lrcorner q -- urcorner q -- ulcorner q --
            reverse  p -- ulcorner q -- cycle
        fi
    fi
    endgroup
enddef ;

boolean intersection_found ;

secondarydef p intersection_point q =
    begingroup
    save x_, y_ ;
    (x_,y_) = p intersectiontimes q ;
    if x_<0 :
        intersection_found := false ;
        center p % origin
    else :
        intersection_found := true ;
        .5[point x_ of p, point y_ of q]
    fi
    endgroup
enddef ;

%D New, undocumented, experimental:

vardef tensecircle (expr width, height, offset) =
    (-width/2,-height/2) ... (0,-height/2-offset) ...
    (+width/2,-height/2) ... (+width/2+offset,0)  ...
    (+width/2,+height/2) ... (0,+height/2+offset) ...
    (-width/2,+height/2) ... (-width/2-offset,0)  ... cycle
enddef ;

vardef roundedsquare (expr width, height, offset) =
    (offset,0)            -- (width-offset,0)      {right} ..
    (width,offset)        -- (width,height-offset) {up}    ..
    (width-offset,height) -- (offset,height)       {left}  ..
    (0,height-offset)     -- (0,offset)            {down}  .. cycle
enddef ;

%D Some colors.

def colortype(expr c) =
    if cmykcolor c : cmykcolor elseif rgbcolor c : rgbcolor else : grayscale fi
enddef ;

vardef whitecolor(expr c) =
    if cmykcolor c : (0,0,0,0) elseif rgbcolor c : (1,1,1)  else : 1         fi
enddef ;

vardef blackcolor(expr c) =
    if cmykcolor c : (0,0,0,1) elseif rgbcolor c : (0,0,0)  else : 0         fi
enddef ;

%D Well, this is the dangerous and naive version:

def drawfill text t =
    fill t ;
    draw t ;
enddef;

%D This two step approach saves the path first, since it can
%D be a function. Attributes must not be randomized.

def drawfill expr c =
    path _c_ ; _c_ := c ;
    mfun_do_drawfill
enddef ;

def mfun_do_drawfill text t =
    draw _c_ t ;
    fill _c_ t ;
enddef;

def undrawfill expr c =
    drawfill c withcolor background % rather useless
enddef ;

%D Moved from mp-char.mp

vardef paired primary d =
    if pair d : d else : (d,d) fi
enddef ;

vardef tripled primary d =
    if color d : d else : (d,d,d) fi
enddef ;

% maybe secondaries:

primarydef p enlarged       d = ( p llmoved d -- p lrmoved d -- p urmoved d -- p ulmoved d -- cycle ) enddef ;
primarydef p llenlarged     d = ( p llmoved d -- lrcorner p  -- urcorner p  -- ulcorner p  -- cycle ) enddef ;
primarydef p lrenlarged     d = ( llcorner p  -- p lrmoved d -- urcorner p  -- ulcorner p  -- cycle ) enddef ;
primarydef p urenlarged     d = ( llcorner p  -- lrcorner p  -- p urmoved d -- ulcorner p  -- cycle ) enddef ;
primarydef p ulenlarged     d = ( llcorner p  -- lrcorner p  -- urcorner p  -- p ulmoved d -- cycle ) enddef ;

primarydef p llmoved        d = ( (llcorner p) shifted (-xpart paired(d),-ypart paired(d)) ) enddef ;
primarydef p lrmoved        d = ( (lrcorner p) shifted (+xpart paired(d),-ypart paired(d)) ) enddef ;
primarydef p urmoved        d = ( (urcorner p) shifted (+xpart paired(d),+ypart paired(d)) ) enddef ;
primarydef p ulmoved        d = ( (ulcorner p) shifted (-xpart paired(d),+ypart paired(d)) ) enddef ;

primarydef p leftenlarged   d = ( (llcorner p) shifted (-d,0) -- lrcorner p -- urcorner p -- (ulcorner p) shifted (-d,0) -- cycle ) enddef ;
primarydef p rightenlarged  d = ( llcorner p -- (lrcorner p) shifted (d,0) -- (urcorner p) shifted (d,0) -- ulcorner p -- cycle   ) enddef ;
primarydef p topenlarged    d = ( llcorner p -- lrcorner p -- (urcorner p) shifted (0,d) -- (ulcorner p) shifted (0,d) -- cycle   ) enddef ;
primarydef p bottomenlarged d = ( llcorner p shifted (0,-d) -- lrcorner p shifted (0,-d) -- urcorner p -- ulcorner p -- cycle     ) enddef ;

%D Handy for testing/debugging:

primarydef p crossed d = (
    if pair p :
        p shifted (-d, 0) -- p --
        p shifted ( 0,-d) -- p --
        p shifted (+d, 0) -- p --
        p shifted ( 0,+d) -- p -- cycle
    else :
        center p shifted (-d, 0) -- llcorner p --
        center p shifted ( 0,-d) -- lrcorner p --
        center p shifted (+d, 0) -- urcorner p --
        center p shifted ( 0,+d) -- ulcorner p -- cycle
    fi
) enddef ;

%D Also handy (math ladders):

vardef laddered primary p = % was expr
    point 0 of p
    for i=1 upto length(p) :
        -- (xpart (point i of p), ypart (point (i-1) of p)) -- (point i of p)
    endfor
enddef ;

%D Saves typing:

% vardef bottomboundary primary p = (llcorner p -- lrcorner p) enddef ;
% vardef rightboundary  primary p = (lrcorner p -- urcorner p) enddef ;
% vardef topboundary    primary p = (urcorner p -- ulcorner p) enddef ;
% vardef leftboundary   primary p = (ulcorner p -- llcorner p) enddef ;

vardef bottomboundary primary p = if pair p : p else : (llcorner p -- lrcorner p) fi enddef ;
vardef rightboundary  primary p = if pair p : p else : (lrcorner p -- urcorner p) fi enddef ;
vardef topboundary    primary p = if pair p : p else : (urcorner p -- ulcorner p) fi enddef ;
vardef leftboundary   primary p = if pair p : p else : (ulcorner p -- llcorner p) fi enddef ;

%D Nice too:

primarydef p superellipsed s =
    superellipse (
        .5[lrcorner p,urcorner p],
        .5[urcorner p,ulcorner p],
        .5[ulcorner p,llcorner p],
        .5[llcorner p,lrcorner p],
        s
    )
enddef ;

primarydef p squeezed s = (
    (llcorner p .. .5[llcorner p,lrcorner p] shifted ( 0, ypart paired(s)) .. lrcorner p) &
    (lrcorner p .. .5[lrcorner p,urcorner p] shifted (-xpart paired(s), 0) .. urcorner p) &
    (urcorner p .. .5[urcorner p,ulcorner p] shifted ( 0,-ypart paired(s)) .. ulcorner p) &
    (ulcorner p .. .5[ulcorner p,llcorner p] shifted ( xpart paired(s), 0) .. llcorner p) & cycle
) enddef ;

primarydef p randomshifted s =
    begingroup ;
    save ss ; pair ss ;
    ss := paired(s) ;
    p shifted (-.5xpart ss + uniformdeviate xpart ss,-.5ypart ss + uniformdeviate ypart ss)
    endgroup
enddef ;

primarydef p randomized s = (
    if path p :
        for i=0 upto length(p)-1 :
            ((point       i    of p) randomshifted s) .. controls
            ((postcontrol i    of p) randomshifted s) and
            ((precontrol (i+1) of p) randomshifted s) ..
        endfor
        if cycle p :
            cycle
        else :
            ((point length(p) of p) randomshifted s)
        fi
    elseif pair p :
        p randomshifted s
    elseif cmykcolor p :
        if color s :
           ((uniformdeviate cyanpart    s) * cyanpart    p,
            (uniformdeviate magentapart s) * magentapart p,
            (uniformdeviate yellowpart  s) * yellowpart  p,
            (uniformdeviate blackpart   s) * blackpart   p)
        elseif pair s :
            ((xpart s + (uniformdeviate (ypart s - xpart s))) * p)
        else :
            ((uniformdeviate s) * p)
        fi
    elseif rgbcolor p :
        if color s :
           ((uniformdeviate redpart   s) * redpart   p,
            (uniformdeviate greenpart s) * greenpart p,
            (uniformdeviate bluepart  s) * bluepart  p)
        elseif pair s :
           ((xpart s + (uniformdeviate (ypart s - xpart s))) * p)
        else :
           ((uniformdeviate s) * p)
        fi
    elseif color p :
        if color s :
            ((uniformdeviate greypart s) * greypart p)
        elseif pair s :
            ((xpart s + (uniformdeviate (ypart s - xpart s))) * p)
        else :
            ((uniformdeviate s) * p)
        fi
    else :
        p + uniformdeviate s
    fi
) enddef ;

%D Not perfect (alternative for interpath)

vardef interpolated(expr s, p, q) =
    save m ; numeric m ;
    m := max(length(p),length(q)) ;
    if path p :
        for i=0 upto m-1 :
            s[point       (i   /m) along p,point       (i   /m) along q] .. controls
            s[postcontrol (i   /m) along p,postcontrol (i   /m) along q] and
            s[precontrol ((i+1)/m) along p,precontrol ((i+1)/m) along q] ..
        endfor
        if cycle p :
            cycle
        else :
            s[point infinity of p,point infinity of q]
        fi
    else :
        a[p,q]
    fi
enddef ;

%D Interesting too:

primarydef p paralleled d = (
    p shifted if d < 0 : - fi ((point abs(d) on (p rotatedaround(point 0 of p,90))) - point 0 of p)
) enddef ;

vardef punked primary p =
    point 0 of p for i=1 upto length(p)-1 : -- point i of p endfor
    if cycle p : -- cycle else : -- point length(p) of p fi
enddef ;

vardef curved primary p =
    point 0 of p for i=1 upto length(p)-1 : .. point i of p endfor
    if cycle p : .. cycle else : .. point length(p) of p fi
enddef ;

primarydef p blownup s =
    begingroup
        save _p_ ; path _p_ ;
        _p_ := p xysized (bbwidth(p)+2(xpart paired(s)),bbheight(p)+2(ypart paired(s))) ;
        (_p_ shifted (center p - center _p_))
    endgroup
enddef ;

%D Rather fundamental.

% not yet ok

vardef leftrightpath(expr p, l) = % used in s-pre-19
    save q, r, t, b ; path q, r ; pair t, b ;
    t := (ulcorner p -- urcorner p) intersection_point p ;
    b := (llcorner p -- lrcorner p) intersection_point p ;
    r := if xpart directionpoint t of p < 0 : reverse p else : p fi ; % r is needed, else problems when reverse is fed
    q := r cutbefore if l: t else: b fi ;
    q := q if xpart point 0 of r > 0 : & r fi cutafter  if l: b else: t fi ;
    q
enddef ;

vardef  leftpath expr p = leftrightpath(p,true ) enddef ;
vardef rightpath expr p = leftrightpath(p,false) enddef ;

%D Drawoptions

def saveoptions =
    save _op_ ; def _op_ = enddef ;
enddef ;

%D Tracing. (not yet in lexer)

let normaldraw = draw ;
let normalfill = fill ;

% bugged in mplib so ...

def normalfill expr c = addto currentpicture contour c _op_ enddef ;
def normaldraw expr p = addto currentpicture if picture p: also p else: doublepath p withpen currentpen fi _op_ enddef ;

def drawlineoptions   (text t) = def _lin_opt_ = t enddef ; enddef ;
def drawpointoptions  (text t) = def _pnt_opt_ = t enddef ; enddef ;
def drawcontroloptions(text t) = def _ctr_opt_ = t enddef ; enddef ;
def drawlabeloptions  (text t) = def _lab_opt_ = t enddef ; enddef ;
def draworiginoptions (text t) = def _ori_opt_ = t enddef ; enddef ;
def drawboundoptions  (text t) = def _bnd_opt_ = t enddef ; enddef ;
def drawpathoptions   (text t) = def _pth_opt_ = t enddef ; enddef ;

def resetdrawoptions =
    drawlineoptions   (withpen pencircle scaled 1pt   withcolor .5white) ;
    drawpointoptions  (withpen pencircle scaled 4pt   withcolor   black) ;
    drawcontroloptions(withpen pencircle scaled 2.5pt withcolor   black) ;
    drawlabeloptions  () ;
    draworiginoptions (withpen pencircle scaled 1pt   withcolor .5white) ;
    drawboundoptions  (dashed evenly _ori_opt_) ;
    drawpathoptions   (withpen pencircle scaled 5pt   withcolor .8white) ;
enddef ;

resetdrawoptions ;

%D Path.

def drawpath expr p =
    normaldraw p _pth_opt_
enddef ;

%D Arrow.

vardef drawarrowpath expr p =
    save autoarrows ; boolean autoarrows ; autoarrows := true ;
    drawarrow p _pth_opt_
enddef ;

% def drawarrowpath expr p =
%   begingroup ;
%   save autoarrows ; boolean autoarrows ; autoarrows := true ;
%   save arrowpath ; path arrowpath ; arrowpath := p ;
%   _drawarrowpath_
% enddef ;
%
% def _drawarrowpath_ text t =
%   drawarrow arrowpath _pth_opt_ t ;
%   endgroup ;
% enddef ;

def midarrowhead expr p =
    arrowhead p cutafter (point length(p cutafter point .5 along p)+ahlength on p)
enddef ;

vardef arrowheadonpath (expr p, s) =
    save autoarrows ; boolean autoarrows ;
    autoarrows := true ;
    set_ahlength(scaled ahfactor) ; % added
    arrowhead p if s<1 : cutafter (point (s*arclength(p)+.5ahlength) on p) fi
enddef ;

%D Points.

def drawpoint expr c =
    if string c :
        string _c_ ;
        _c_ := "(" & c & ")" ;
        dotlabel.urt(_c_, scantokens _c_) ;
        drawdot scantokens _c_
    else :
        dotlabel.urt("(" & decimal xpart c & "," & decimal ypart c & ")", c) ;
        drawdot c
    fi _pnt_opt_
enddef ;

%D PathPoints.

def drawpoints        expr c = path _c_ ; _c_ := c ; mfun_draw_points        enddef ;
def drawcontrolpoints expr c = path _c_ ; _c_ := c ; mfun_draw_controlpoints enddef ;
def drawcontrollines  expr c = path _c_ ; _c_ := c ; mfun_draw_controllines  enddef ;
def drawpointlabels   expr c = path _c_ ; _c_ := c ; mfun_draw_pointlabels   enddef ;

def mfun_draw_points text t =
    for _i_=0 upto length(_c_) :
        normaldraw point _i_ of _c_ _pnt_opt_ t ;
    endfor ;
enddef;

def mfun_draw_controlpoints text t =
    for _i_=0 upto length(_c_) :
        normaldraw precontrol  _i_ of _c_ _ctr_opt_ t ;
        normaldraw postcontrol _i_ of _c_ _ctr_opt_ t ;
    endfor ;
enddef;

def mfun_draw_controllines text t =
    for _i_=0 upto length(_c_) :
        normaldraw point _i_ of _c_ -- precontrol  _i_ of _c_ _lin_opt_ t ;
        normaldraw point _i_ of _c_ -- postcontrol _i_ of _c_ _lin_opt_ t ;
    endfor ;
enddef;

boolean swappointlabels ; swappointlabels := false ;

def mfun_draw_pointlabels text t =
    for _i_=0 upto length(_c_) :
        pair _u_ ; _u_ := unitvector(direction _i_ of _c_) rotated if swappointlabels : - fi 90 ;
        pair _p_ ; _p_ := (point _i_ of _c_) ;
        _u_ := 12 * defaultscale * _u_ ;
        normaldraw thelabel ( decimal _i_, _p_ shifted if cycle _c_ and (_i_=0) : - fi _u_ ) _lab_opt_ t ;
    endfor ;
enddef;

%D Bounding box.

def drawboundingbox expr p =
    normaldraw boundingbox p _bnd_opt_
enddef ;

%D Origin.

numeric originlength ; originlength := .5cm ;

def draworigin text t =
    normaldraw (origin shifted (0, originlength) -- origin shifted (0,-originlength)) _ori_opt_ t ;
    normaldraw (origin shifted ( originlength,0) -- origin shifted (-originlength,0)) _ori_opt_ t ;
enddef;

%D Axis.

numeric tickstep   ; tickstep   := 5mm ;
numeric ticklength ; ticklength := 2mm ;

def drawxticks expr c = path _c_ ; _c_ := c ; mfun_draw_xticks enddef ;
def drawyticks expr c = path _c_ ; _c_ := c ; mfun_draw_yticks enddef ;
def drawticks  expr c = path _c_ ; _c_ := c ; mfun_draw_ticks  enddef ;

% Adding eps prevents disappearance due to rounding errors.

def mfun_draw_xticks text t =
    for i=0 step -tickstep until xpart llcorner _c_ - eps :
        if (i<=xpart lrcorner _c_) :
        normaldraw (i,-ticklength)--(i,ticklength) _ori_opt_ t ;
        fi ;
    endfor ;
    for i=0 step  tickstep until xpart lrcorner _c_ + eps :
        if (i>=xpart llcorner _c_) :
            normaldraw (i,-ticklength)--(i,ticklength) _ori_opt_ t ;
        fi ;
    endfor ;
    normaldraw (llcorner _c_ -- ulcorner _c_) shifted (-xpart llcorner _c_,0) _ori_opt_ t ;
enddef ;

def mfun_draw_yticks text t =
    for i=0 step -tickstep until ypart llcorner _c_ - eps :
        if (i<=ypart ulcorner _c_) :
            normaldraw (-ticklength,i)--(ticklength,i) _ori_opt_ t ;
        fi ;
    endfor ;
    for i=0 step  tickstep until ypart ulcorner _c_ + eps :
        if (i>=ypart llcorner _c_) :
            normaldraw (-ticklength,i)--(ticklength,i) _ori_opt_ t ;
        fi ;
    endfor ;
    normaldraw (llcorner _c_ -- lrcorner _c_) shifted (0,-ypart llcorner _c_) _ori_opt_ t ;
enddef ;

def mfun_draw_ticks text t =
    drawxticks _c_ t ;
    drawyticks _c_ t ;
enddef ;

%D All of it except axis.

def drawwholepath expr p =
    draworigin          ;
    drawpath          p ;
    drawcontrollines  p ;
    drawcontrolpoints p ;
    drawpoints        p ;
    drawboundingbox   p ;
    drawpointlabels   p ;
enddef ;

%D Tracing.

def visualizeddraw expr c =
    if picture c : normaldraw c else : path _c_ ; _c_ := c ; do_visualizeddraw fi
enddef ;

def visualizedfill expr c =
    if picture c : normalfill c else : path _c_ ; _c_ := c ; do_visualizedfill fi
enddef ;

def do_visualizeddraw text t =
    draworigin            ;
    drawpath          _c_ t ;
    drawcontrollines  _c_ ;
    drawcontrolpoints _c_ ;
    drawpoints        _c_ ;
    drawboundingbox   _c_ ;
    drawpointlabels   _c_ ;
enddef ;

def do_visualizedfill text t =
    if cycle _c_ : normalfill _c_ t fi ;
    draworigin            ;
    drawcontrollines  _c_ ;
    drawcontrolpoints _c_ ;
    drawpoints        _c_ ;
    drawboundingbox   _c_ ;
    drawpointlabels   _c_ ;
enddef ;

def visualizepaths =
    let fill = visualizedfill ;
    let draw = visualizeddraw ;
enddef ;

def naturalizepaths =
    let fill = normalfill ;
    let draw = normaldraw ;
enddef ;

extra_endfig := extra_endfig & " naturalizepaths ; " ;

%D Nice tracer:

def drawboundary primary p =
    draw p dashed evenly withcolor white ;
    draw p dashed oddly  withcolor black ;
    draw (- llcorner p) withpen pencircle scaled 3   withcolor white ;
    draw (- llcorner p) withpen pencircle scaled 1.5 withcolor black ;
enddef ;

%D Also handy:

extra_beginfig := extra_beginfig & " truecorners :=  0 ; "   ; % restores
extra_beginfig := extra_beginfig & " miterlimit  := 10 ; "   ; % restores
extra_beginfig := extra_beginfig & " linejoin := rounded ; " ; % restores
extra_beginfig := extra_beginfig & " linecap  := rounded ; " ; % restores

%D Normally, arrowheads don't scale well. So we provide a
%D hack.

boolean autoarrows ; autoarrows := false ;
numeric ahfactor   ; ahfactor   := 2.5 ;

def set_ahlength (text t) =
  % ahlength := (ahfactor*pen_size(_op_ t)) ; % _op_ added
  % problem: _op_ can contain color so a no-go, we could apply the transform
  % but i need to figure out the best way (fakepicture and take components).
    ahlength := (ahfactor*pen_size(t)) ;
enddef ;

vardef pen_size (text t) =
    save p ; picture p ; p := nullpicture ;
    addto p doublepath (top origin -- bot origin) t ;
    (ypart urcorner p - ypart lrcorner p)
enddef ;

%D The next two macros are adapted versions of plain
%D \METAPOST\ definitions.

vardef arrowpath expr p = % patch by Peter Rolf: supports squared pen and shifting (hh: maybe just use center of head as first)
    (p cutafter makepath(pencircle scaled 2(ahlength*cosd(.5ahangle)) shifted point length p of p))
enddef;

% def _finarr text t =
%     if autoarrows : set_ahlength (t) fi ;
%     draw     arrowpath _apth t ; % arrowpath added
%     filldraw arrowhead _apth t ;
% enddef;

def _finarr text t =
    if autoarrows : set_ahlength (t) fi ;
    draw arrowpath _apth t ; % arrowpath added
    fill arrowhead _apth t ;
    draw arrowhead _apth t ;
enddef;

def _finarr text t =
    if autoarrows : set_ahlength (t) fi ;
    draw arrowpath _apth t ; % arrowpath added
    fill arrowhead _apth t ;
    draw arrowhead _apth t undashed ;
enddef;

%D Handy too ......

vardef pointarrow (expr pat, loc, len, off) =
    save l, r, s, t ; path l, r ; numeric s ; pair t ;
    t := if pair loc : loc else : point loc along pat fi ;
    s := len/2 - off ; if s<=0 : s := 0 elseif s>len : s := len fi ;
    r := pat cutbefore t ;
    r := (r cutafter point (arctime s of r) of r) ;
    s := len/2 + off ; if s<=0 : s := 0 elseif s>len : s := len fi ;
    l := reverse (pat cutafter t) ;
    l := (reverse (l cutafter point (arctime s of l) of l)) ;
    (l..r)
enddef ;

def rightarrow  (expr pat,tim,len) = pointarrow(pat,tim,len,-len) enddef ;
def leftarrow   (expr pat,tim,len) = pointarrow(pat,tim,len,+len) enddef ;
def centerarrow (expr pat,tim,len) = pointarrow(pat,tim,len,   0) enddef ;

%D The \type {along} and \type {on} operators can be used
%D as follows:
%D
%D \starttyping
%D drawdot point .5  along somepath ;
%D drawdot point 3cm on    somepath ;
%D \stoptyping
%D
%D The number denotes a percentage (fraction).

primarydef pct along pat = % also negative
    (arctime (pct * (arclength pat)) of pat) of pat
enddef ;

primarydef len on pat = % no outer ( ) .. somehow fails
    (arctime if len>=0 : len else : (arclength(pat)+len) fi of pat) of pat
enddef ;

% this cuts of a piece from both ends

% tertiarydef pat cutends len =
%   begingroup ; save tap ; path tap ;
%   tap := pat cutbefore (point len on pat) ;
%   (tap cutafter (point -len on tap))
%   endgroup
% enddef ;

tertiarydef pat cutends len =
    begingroup
    save tap ; path tap ;
    tap := pat cutbefore (point (xpart paired(len)) on pat) ;
    (tap cutafter (point -(ypart paired(len)) on tap))
    endgroup
enddef ;

%D To be documented.

path freesquare ;

freesquare := (
    (-1,0) -- (-1,-1) -- (0,-1) -- (+1,-1) --
    (+1,0) -- (+1,+1) -- (0,+1) -- (-1,+1) -- cycle
) scaled .5 ;

numeric freelabeloffset  ; freelabeloffset  := 3pt ;
numeric freedotlabelsize ; freedotlabelsize := 3pt ;

vardef thefreelabel (expr str, loc, ori) =
    save s, p, q, l ; picture s ; path p, q ; pair l ;
    interim labeloffset := freelabeloffset ;
    s := if string str : thelabel(str,loc) else : str shifted -center str shifted loc fi ;
    setbounds s to boundingbox s enlarged freelabeloffset ;
    p := fullcircle scaled (2*length(loc-ori)) shifted ori ;
    q := freesquare xyscaled (urcorner s - llcorner s) ;
    l := point xpart (p intersectiontimes (ori--loc shifted (loc-ori))) of q ;
    setbounds s to boundingbox s enlarged -freelabeloffset ; % new
  % draw boundingbox s shifted -l withpen pencircle scaled .5pt withcolor red ;
    (s shifted -l)
enddef ;

vardef freelabel (expr str, loc, ori) =
    draw thefreelabel(str,loc,ori) ;
enddef ;

vardef freedotlabel (expr str, loc, ori) =
    interim linecap := rounded ;
    draw loc withpen pencircle scaled freedotlabelsize ;
    draw thefreelabel(str,loc,ori) ;
enddef ;

%D \starttyping
%D drawarrow anglebetween(line_a,line_b,somelabel) ;
%D \stoptyping

newinternal angleoffset ; angleoffset :=  0pt ;
newinternal anglelength ; anglelength := 20pt ;
newinternal anglemethod ; anglemethod :=    1 ;

% vardef anglebetween (expr a, b, str) = % path path string
%   save pointa, pointb, common, middle, offset ;
%   pair pointa, pointb, common, middle, offset ;
%   save curve ; path curve ;
%   save where ; numeric where ;
%   if round point 0 of a = round point 0 of b :
%     common := point 0 of a ;
%   else :
%     common := a intersectionpoint b ;
%   fi ;
%   pointa := point anglelength on a ;
%   pointb := point anglelength on b ;
%   where  := turningnumber (common--pointa--pointb--cycle) ;
%   middle := ((common--pointa) rotatedaround (pointa,-where*90))
%                             intersectionpoint
%             ((common--pointb) rotatedaround (pointb, where*90)) ;
%   if     anglemethod = 0 :
%     curve  := pointa{unitvector(middle-pointa)}.. pointb;
%     middle := point .5 along curve ;
%     curve  := common ;
%   elseif anglemethod = 1 :
%     curve  := pointa{unitvector(middle-pointa)}.. pointb;
%     middle := point .5 along curve ;
%   elseif anglemethod = 2 :
%     middle := common rotatedaround(.5[pointa,pointb],180) ;
%     curve  := pointa--middle--pointb ;
%   elseif anglemethod = 3 :
%     curve  := pointa--middle--pointb ;
%   elseif anglemethod = 4 :
%     curve  := pointa..controls middle..pointb ;
%     middle := point .5 along curve ;
%   fi ;
%   draw thefreelabel(str, middle, common) withcolor black ;
%   curve
% enddef ;

vardef anglebetween (expr a, b, str) = % path path string
    save pointa, pointb, common, middle, offset ;
    pair pointa, pointb, common, middle, offset ;
    save curve ; path curve ;
    save where ; numeric where ;
    if round point 0 of a = round point 0 of b :
        common := point 0 of a ;
    else :
        common := a intersectionpoint b ;
    fi ;
    pointa := point anglelength on a ;
    pointb := point anglelength on b ;
    where  := turningnumber (common--pointa--pointb--cycle) ;
    middle := (reverse(common--pointa) rotatedaround (pointa,-where*90))
        intersection_point
        (reverse(common--pointb) rotatedaround (pointb, where*90)) ;
    if not intersection_found :
        middle := point .5 along
        ((reverse(common--pointa) rotatedaround (pointa,-where*90)) --
        (       (common--pointb) rotatedaround (pointb, where*90))) ;
    fi ;
    if     anglemethod = 0 :
        curve  := pointa{unitvector(middle-pointa)}.. pointb;
        middle := point .5 along curve ;
        curve  := common ;
    elseif anglemethod = 1 :
        curve  := pointa{unitvector(middle-pointa)}.. pointb;
        middle := point .5 along curve ;
    elseif anglemethod = 2 :
        middle := common rotatedaround(.5[pointa,pointb],180) ;
        curve  := pointa--middle--pointb ;
    elseif anglemethod = 3 :
        curve  := pointa--middle--pointb ;
    elseif anglemethod = 4 :
        curve  := pointa..controls middle..pointb ;
        middle := point .5 along curve ;
    fi ;
    draw thefreelabel(str, middle, common) ; % withcolor black ;
    curve
enddef ;

% Stack

picture mfun_current_picture_stack[] ;
numeric mfun_current_picture_depth   ;

mfun_current_picture_depth := 0 ;

def pushcurrentpicture =
    mfun_current_picture_depth := mfun_current_picture_depth + 1 ;
    mfun_current_picture_stack[mfun_current_picture_depth] := currentpicture ;
    currentpicture := nullpicture ;
enddef ;

def popcurrentpicture text t = % optional text
    if mfun_current_picture_depth > 0 :
        addto mfun_current_picture_stack[mfun_current_picture_depth] also currentpicture t ;
        currentpicture := mfun_current_picture_stack[mfun_current_picture_depth] ;
        mfun_current_picture_stack[mfun_current_picture_depth] := nullpicture ;
        mfun_current_picture_depth := mfun_current_picture_depth - 1 ;
    fi ;
enddef ;

%D colorcircle(size, red, green, blue) ;

% vardef colorcircle (expr size, red, green, blue) =
%   save r, g, b, rr, gg, bb, cc, mm, yy ; save radius ;
%   path r, g, b, rr, bb, gg, cc, mm, yy ; numeric radius ;
%
%   radius := 5cm ; pickup pencircle scaled (radius/25) ;
%
%   r := g := b := fullcircle scaled radius shifted (0,radius/4) ;
%
%   r := r rotatedaround (origin, 15) ;
%   g := g rotatedaround (origin,135) ;
%   b := b rotatedaround (origin,255) ;
%
%   r := r rotatedaround(center r,-90) ;
%   g := g rotatedaround(center g, 90) ;
%
%   gg := buildcycle(buildcycle(reverse r,b),g) ;
%   cc := buildcycle(buildcycle(b,reverse g),r) ;
%
%   rr := gg rotatedaround(origin,120) ;
%   bb := gg rotatedaround(origin,240) ;
%
%   yy := cc rotatedaround(origin,120) ;
%   mm := cc rotatedaround(origin,240) ;
%
%   pushcurrentpicture ;
%
%   fill fullcircle scaled radius withcolor white ;
%
%   fill rr withcolor red   ; fill cc withcolor white-red   ;
%   fill gg withcolor green ; fill mm withcolor white-green ;
%   fill bb withcolor blue  ; fill yy withcolor white-blue  ;
%
%   for i = rr,gg,bb,cc,mm,yy : draw i withcolor .5white ; endfor ;
%
%   currentpicture := currentpicture xsized size ;
%
%   popcurrentpicture ;
% enddef ;

% vardef colorcircle (expr size, red, green, blue) =
%   save r, g, b, rr, gg, bb, cc, mm, yy ; save radius ;
%   path r, g, b, rr, bb, gg, cc, mm, yy ; numeric radius ;
%
%   radius := 5cm ; pickup pencircle scaled (radius/25) ;
%
%   transform t ; t := identity rotatedaround(origin,120) ;
%
%   r := fullcircle scaled radius
%    shifted (0,radius/4) rotatedaround(origin,15) ;
%
%   g := r transformed t ; b := g transformed t ;
%
%   r := r rotatedaround(center r,-90) ;
%   g := g rotatedaround(center g, 90) ;
%
%   gg := buildcycle(buildcycle(reverse r,b),g) ;
%   cc := buildcycle(buildcycle(b,reverse g),r) ;
%
%   rr := gg transformed t ; bb := rr transformed t ;
%   yy := cc transformed t ; mm := yy transformed t ;
%
%   pushcurrentpicture ;
%
%   fill fullcircle scaled radius withcolor white ;
%
%   fill rr withcolor red   ; fill cc withcolor white-red   ;
%   fill gg withcolor green ; fill mm withcolor white-green ;
%   fill bb withcolor blue  ; fill yy withcolor white-blue  ;
%
%   for i = rr,gg,bb,cc,mm,yy : draw i withcolor .5white ; endfor ;
%
%   currentpicture := currentpicture xsized size ;
%
%   popcurrentpicture ;
% enddef ;

vardef colorcircle (expr size, red, green, blue) = % might move
    save r, g, b, c, m, y, w ; save radius ;
    path r, g, b, c, m, y, w ; numeric radius ;

    radius := 5cm ; pickup pencircle scaled (radius/25) ;

    transform t ; t := identity rotatedaround(origin,120) ;

    r := fullcircle rotated 90 scaled radius shifted (0,radius/4) rotatedaround(origin,135) ;

    b := r transformed t ; g := b transformed t ;

    c := buildcycle(subpath(1,7) of g,subpath(1,7) of b) ;
    y := c transformed t ; m := y transformed t ;

    w := buildcycle(subpath(3,5) of r, subpath(3,5) of g,subpath(3,5) of b) ;

    pushcurrentpicture ;

    fill r withcolor         red   ;
    fill g withcolor         green ;
    fill b withcolor         blue  ;
    fill c withcolor white - red   ;
    fill m withcolor white - green ;
    fill y withcolor white - blue  ;
    fill w withcolor white         ;

    for i = r,g,b,c,m,y : draw i withcolor .5white ; endfor ;

    currentpicture := currentpicture xsized size ;

    popcurrentpicture ;
enddef ;

% penpoint (i,2) of somepath -> inner / outer point

vardef penpoint expr pnt of p =
    save n, d ; numeric n, d ;
    (n,d) = if pair pnt : pnt else : (pnt,1) fi ;
    (point n of p shifted ((penoffset direction n of p of currentpen) scaled d))
enddef ;

% nice: currentpicture := inverted currentpicture ;

primarydef p uncolored c =
    if color p :
        c - p
    else :
        image (
            for i within p :
                addto currentpicture
                    if stroked i or filled i :
                        if filled i :
                            contour
                        else :
                            doublepath
                        fi
                        pathpart i
                        dashed dashpart i withpen penpart i
                    else :
                        also i
                    fi
                    withcolor c-(redpart i, greenpart i, bluepart i) ;
            endfor ;
        )
  fi
enddef ;

vardef inverted primary p =
    p uncolored white
enddef ;

% primarydef p softened c =
%   if color p :
%     tripled(c) * p
%   else :
%     image
%       (save cc ; color cc ; cc := tripled(c) ;
%        for i within p :
%          addto currentpicture
%          if stroked i or filled i :
%            if filled i : contour else : doublepath fi pathpart i
%            dashed dashpart i withpen penpart i
%          else :
%           also i
%          fi
%          withcolor (redpart   cc * redpart   i,
%                     greenpart cc * greenpart i,
%                     bluepart  cc * bluepart  i) ;
%        endfor ;)
%   fi
% enddef ;

primarydef p softened c =
    begingroup
    save cc ; color cc ; cc := tripled(c) ;
    if color p :
        (redpart cc * redpart p,greenpart cc * greenpart p, bluepart cc * bluepart p)
    else :
        image (
            for i within p :
                addto currentpicture
                    if stroked i or filled i :
                        if filled i :
                            contour
                        else :
                            doublepath
                        fi
                        pathpart i
                        dashed dashpart i withpen penpart i
                    else :
                        also i
                    fi
                    withcolor (redpart cc * redpart i, greenpart cc * greenpart i, bluepart cc * bluepart i) ;
            endfor ;
        )
    fi
    endgroup
enddef ;

vardef grayed primary p =
    if color p :
        tripled(.30redpart p+.59greenpart p+.11bluepart p)
    else :
        image (
            for i within p :
                addto currentpicture
                    if stroked i or filled i :
                        if filled i :
                            contour
                        else :
                            doublepath
                        fi
                        pathpart i
                        dashed dashpart i
                        withpen penpart i
                    else :
                        also i
                    fi
                withcolor tripled(.30redpart i+.59greenpart i+.11bluepart i) ;
            endfor ;
        )
  fi
enddef ;

% yes or no: "text" infont "cmr12" at 24pt ;

% let normalinfont = infont ;
%
% numeric lastfontsize ; lastfontsize = fontsize defaultfont ;
%
% def infont primary name =  % no vardef, no expr
%   hide(lastfontsize := fontsize name) % no ;
%   normalinfont name
% enddef ;
%
% def scaledat expr size =
%   scaled (size/lastfontsize)
% enddef ;
%
% let at = scaledat ;

% like decimal

def condition primary b = if b : "true" else : "false" fi enddef ;

% undocumented

primarydef p stretched s =
    begingroup
    save pp ; path pp ; pp := p xyscaled s ;
    (pp shifted ((point 0 of p) - (point 0 of pp)))
    endgroup
enddef ;

% primarydef p enlonged len =
%     begingroup
%     save al ; al := arclength(p) ;
%     if al > 0 :
%         if pair p :
%             point 1 of ((origin -- p) stretched ((al+len)/al))
%         else :
%             p stretched ((al+len)/al)
%         fi
%     else :
%         p
%     fi
%     endgroup
% enddef ;

primarydef p enlonged len =
    begingroup
    if len == 0 :
        p
    elseif pair p :
        save q ; path q ; q := origin -- p ;
        save al ; al := arclength(q) ;
        if al > 0 :
            point 1 of (q stretched ((al+len)/al))
        else :
            p
        fi
    else :
        save al ; al := arclength(p) ;
        if al > 0 :
            p stretched ((al+len)/al)
        else :
            p
        fi
    fi
    endgroup
enddef ;

% path p ; p := (0,0) -- (10cm,5cm) ;
% drawarrow p withcolor red ;
% drawarrow p shortened 1cm withcolor green ;

% primarydef p shortened d =
%     reverse ( ( reverse (p enlonged -d) ) enlonged -d )
% enddef ;

primarydef p shortened d =
    reverse ( ( reverse (p enlonged -xpart paired(d)) ) enlonged -ypart paired(d) )
enddef ;

% yes or no, untested -)

def xshifted expr dx = shifted(dx,0) enddef ;
def yshifted expr dy = shifted(0,dy) enddef ;

% also handy

% right: str = readfrom ("abc" & ".def" ) ;
% wrong: str = readfrom  "abc" & ".def"   ;

% Every 62th read fails so we need to try again!

% def readfile (expr name) =
%   if (readfrom (name) <> EOF) :
%     scantokens("input " & name & ";") ;
%   elseif (readfrom (name) <> EOF) :
%     scantokens("input " & name & ";") ;
%   fi ;
%   closefrom (name) ;
% enddef ;
%
% this sometimes fails on the elseif, so :
%

def readfile (expr name) =
    begingroup ; save ok ; boolean ok ;
    if (readfrom (name) <> EOF) :
        ok := false ;
    elseif (readfrom (name) <> EOF) :
        ok := false ;
    else :
        ok := true ;
    fi ;
    if not ok :
        scantokens("input " & name & " ") ;
    fi ;
    closefrom (name) ;
    endgroup ;
enddef ;

% permits redefinition of end in macro

inner end ;

% this will be redone (when needed) using scripts and backend handling

let normalwithcolor = withcolor ;

def remapcolors =
    def withcolor primary c = normalwithcolor remappedcolor(c) enddef ;
enddef ;

def normalcolors =
    let withcolor = normalwithcolor ;
enddef ;

def resetcolormap =
    color color_map[][][] ;
    normalcolors ;
enddef ;

resetcolormap ;

% color_map_resolution := 1000 ;
%
% def r_color primary c = round(color_map_resolution*redpart   c) enddef ;
% def g_color primary c = round(color_map_resolution*greenpart c) enddef ;
% def b_color primary c = round(color_map_resolution*bluepart  c) enddef ;

def r_color primary c = redpart   c enddef ;
def g_color primary c = greenpart c enddef ;
def b_color primary c = bluepart  c enddef ;

def remapcolor(expr old, new) =
    color_map[redpart old][greenpart old][bluepart old] := new ;
enddef ;

def remappedcolor(expr c) =
    if known color_map[redpart c][greenpart c][bluepart c] :
        color_map[redpart c][greenpart c][bluepart c]
    else :
        c
    fi
enddef ;

% def refill  suffix c = do_repath (1) (c) enddef ;
% def redraw  suffix c = do_repath (2) (c) enddef ;
% def recolor suffix c = do_repath (0) (c) enddef ;
%
% color refillbackground ; refillbackground := (1,1,1) ;
%
% def do_repath (expr mode) (suffix c) text t = % can it be stroked and filled at the same time ?
%   begingroup ;
%   if mode=0 : save withcolor ; remapcolors ; fi ;
%   save _c_, _cc_, _f_, _b_ ; picture _c_, _cc_ ; color _f_ ; path _b_ ;
%   _c_ := c ; _b_ := boundingbox c ; c := nullpicture ;
%   for i within _c_ :
%     _f_ := (redpart i, greenpart i, bluepart i) ;
%     if     bounded i :
%       setbounds c to pathpart i ;
%     elseif clipped i :
%       clip c to pathpart i ;
%     elseif stroked i :
%       addto c doublepath pathpart i
%         dashed dashpart i withpen penpart i
%         withcolor _f_ % (redpart i, greenpart i, bluepart i)
%         if mode=2 : t fi ;
%     elseif filled  i :
%       addto c contour pathpart i
%         withcolor _f_
%         if (mode=1) and (_f_<>refillbackground) : t fi ;
%     else :
%       addto c also i ;
%     fi ;
%   endfor ;
%   setbounds c to _b_ ;
%   endgroup ;
% enddef ;

% Thanks to Jens-Uwe Morawski for pointing out that we need
% to treat bounded and clipped components as local pictures.

def recolor suffix p = p := repathed (0,p) enddef ;
def refill  suffix p = p := repathed (1,p) enddef ;
def redraw  suffix p = p := repathed (2,p) enddef ;
def retext  suffix p = p := repathed (3,p) enddef ;
def untext  suffix p = p := repathed (4,p) enddef ;

% primarydef p recolored t = repathed(0,p) t enddef ;
% primarydef p refilled  t = repathed(1,p) t enddef ;
% primarydef p redrawn   t = repathed(2,p) t enddef ;
% primarydef p retexted  t = repathed(3,p) t enddef ;
% primarydef p untexted  t = repathed(4,p) t enddef ;

color refillbackground ; refillbackground := (1,1,1) ;

% vardef repathed (expr mode, p) text t =
%   begingroup ;
%   if mode=0 : save withcolor ; remapcolors ; fi ;
%   save _p_, _pp_, _f_, _b_, _t_ ;
%   picture _p_, _pp_ ; color _f_ ; path _b_ ; transform _t_ ;
%   _b_ := boundingbox p ; _p_ := nullpicture ;
%   for i within p :
%     _f_ := (redpart i, greenpart i, bluepart i) ;
%     if     bounded i :
%       _pp_ := repathed(mode,i) t ;
%       setbounds _pp_ to pathpart i ;
%       addto _p_ also _pp_ ;
%     elseif clipped i :
%       _pp_ := repathed(mode,i) t ;
%       clip _pp_ to pathpart i ;
%       addto _p_ also _pp_ ;
%     elseif stroked i :
%       addto _p_ doublepath pathpart i
%         dashed dashpart i withpen penpart i
%         withcolor _f_ % (redpart i, greenpart i, bluepart i)
%         if mode=2 : t fi ;
%     elseif filled  i :
%       addto _p_ contour pathpart i
%         withcolor _f_
%         if (mode=1) and (_f_<>refillbackground) : t fi ;
%     elseif textual i : % textpart i <> "" :
%       if mode <> 4 :
%         % transform _t_ ;
%         % (xpart _t_, xxpart _t_, xypart _t_)  = (xpart  i, xxpart i, xypart i) ;
%         % (ypart _t_, yypart _t_, yxpart _t_)  = (ypart  i, yypart i, yxpart i) ;
%         % addto _p_ also
%         %   textpart i infont fontpart i % todo : other font
%         %   transformed _t_
%         %   withpen penpart i
%         %   withcolor _f_
%         %   if mode=3 : t fi ;
%         addto _p_ also i if mode=3 : t fi ;
%       fi ;
%     else :
%       addto _p_ also i ;
%     fi ;
%   endfor ;
%   setbounds _p_ to _b_ ;
%   _p_
%   endgroup
% enddef ;

def restroke  suffix p = p := repathed (21,p) enddef ; % keep attributes
def reprocess suffix p = p := repathed (22,p) enddef ; % no attributes

% also 11 and 12

vardef repathed (expr mode, p) text t =
    begingroup ;
    if mode = 0 :
        save withcolor ;
        remapcolors ;
    fi ;
    save _p_, _pp_, _ppp_, _f_, _b_, _t_ ;
    picture _p_, _pp_, _ppp_ ; color _f_ ; path _b_ ; transform _t_ ;
    _b_ := boundingbox p ;
    _p_ := nullpicture ;
    for i within p :
        _f_ := (redpart i, greenpart i, bluepart i) ;
        if bounded i :
            _pp_ := repathed(mode,i) t ;
            setbounds _pp_ to pathpart i ;
            addto _p_ also _pp_ ;
        elseif clipped i :
            _pp_ := repathed(mode,i) t ;
            clip _pp_ to pathpart i ;
            addto _p_ also _pp_ ;
        elseif stroked i :
            if mode=21 :
                _ppp_ := i ; % indirectness is needed
                addto _p_ also image(scantokens(t & " pathpart _ppp_")
                    dashed dashpart i withpen penpart i
                    withcolor _f_ ; ) ;
            elseif mode=22 :
                _ppp_ := i ; % indirectness is needed
                addto _p_ also image(scantokens(t & " pathpart _ppp_")) ;
            else :
                addto _p_ doublepath pathpart i
                    dashed dashpart i withpen penpart i
                    withcolor _f_ % (redpart i, greenpart i, bluepart i)
                    if mode = 2 :
                        t
                    fi ;
            fi ;
        elseif filled  i :
            if mode=11 :
                _ppp_ := i ; % indirectness is needed
                addto _p_ also image(scantokens(t & " pathpart _ppp_")
                    withcolor _f_ ; ) ;
            elseif mode=12 :
                _ppp_ := i ; % indirectness is needed
                addto _p_ also image(scantokens(t & " pathpart _ppp_")) ;
            else :
                addto _p_ contour pathpart i
                    withcolor _f_
                    if (mode=1) and (_f_<>refillbackground) :
                        t
                    fi ;
            fi ;
        elseif textual i : % textpart i <> "" :
            if mode <> 4 :
                % transform _t_ ;
                % (xpart _t_, xxpart _t_, xypart _t_)  = (xpart  i, xxpart i, xypart i) ;
                % (ypart _t_, yypart _t_, yxpart _t_)  = (ypart  i, yypart i, yxpart i) ;
                % addto _p_ also
                %     textpart i infont fontpart i % todo : other font
                %     transformed _t_
                %     withpen penpart i
                %     withcolor _f_
                %     if mode=3 : t fi ;
                addto _p_ also i
                    if mode=3 :
                        t
                    fi ;
            fi ;
        else :
            addto _p_ also i ;
        fi ;
    endfor ;
    setbounds _p_ to _b_ ;
    _p_
    endgroup
enddef ;

% After a question of Denis on how to erase a z variable, Jacko
% suggested to assign whatever to x and y. So a clearz
% variable can be defined as:
%
% vardef clearz@# =
%   x@# := whatever ;
%   y@# := whatever ;
% enddef ;
%
% but Jacko suggested a redefinition of clearxy:
%
% def clearxy text s =
%  clearxy_index_:=0;
%  for $:=s:
%    clearxy_index_:=clearxy_index_+1; endfor;
%  if clearxy_index_=0:
%    save x,y;
%  else:
%    forsuffixes $:=s: x$:=whatever; y$:=whatever; endfor;
%  fi
% enddef;
%
% which i decided to simplify to:

def clearxy text s =
    if false for $ := s : or true endfor :
        forsuffixes $ := s : x$ := whatever ; y$ := whatever ; endfor ;
    else :
        save x, y ;
    fi
enddef ;

% so now we can say: clearxy ; as well as clearxy 1, 2, 3 ;

% show x0 ; z0 = (10,10) ;
% show x0 ; x0 := whatever ; y0 := whatever ;
% show x0 ; z0 = (20,20) ;
% show x0 ; clearxy 0 ;
% show x0 ; z0 = (30,30) ;

primarydef p smoothed d =
   (p llmoved (-xpart paired(d),0) -- p lrmoved (-xpart paired(d),0) {right} ..
    p lrmoved (0,-ypart paired(d)) -- p urmoved (0,-ypart paired(d)) {up}    ..
    p urmoved (-xpart paired(d),0) -- p ulmoved (-xpart paired(d),0) {left}  ..
    p ulmoved (0,-ypart paired(d)) -- p llmoved (0,-ypart paired(d)) {down}  .. cycle)
enddef ;

primarydef p cornered c =
    ((point 0 of p) shifted (c*(unitvector(point 1 of p - point 0 of p))) --
    for i=1 upto length(p) :
        (point i-1 of p) shifted (c*(unitvector(point i   of p - point i-1 of p))) --
        (point i   of p) shifted (c*(unitvector(point i-1 of p - point i   of p))) ..
        controls point i of p ..
    endfor cycle)
enddef ;

% cmyk color support

vardef cmyk(expr c,m,y,k) =
    (1-c-k,1-m-k,1-y-k)
enddef ;

% handy

% vardef bbwidth (expr p) = % vardef width_of primary p =
%     if known p :
%         if path p or picture p :
%             xpart (lrcorner p - llcorner p)
%         else :
%             0
%         fi
%     else :
%         0
%     fi
% enddef ;

vardef bbwidth primary p =
    if unknown p :
        0
    elseif path p or picture p :
        xpart (lrcorner p - llcorner p)
    else :
        0
    fi
enddef ;

% vardef bbheight (expr p) = % vardef heigth_of primary p =
%     if known p :
%         if path p or picture p :
%             ypart (urcorner p - lrcorner p)
%         else :
%             0
%         fi
%     else :
%         0
%     fi
% enddef ;

vardef bbheight primary p =
    if unknown p :
        0
    elseif path p or picture p :
        ypart (urcorner p - lrcorner p)
    else :
        0
    fi
enddef ;

color nocolor ; numeric noline ; % both unknown signals

def dowithpath (expr p, lw, lc, bc) =
    if known p :
        if known bc :
            fill p withcolor bc ;
        fi ;
        if known lw and known lc :
            draw p withpen pencircle scaled lw withcolor lc ;
        elseif known lw :
            draw p withpen pencircle scaled lw ;
        elseif known lc :
            draw p withcolor lc ;
        fi ;
    fi ;
enddef ;

% result from metafont discussion list (denisr/boguslawj)

def [[ = [ [ enddef ; def [[[ = [ [ [ enddef ;
def ]] = ] ] enddef ; def ]]] = ] ] ] enddef ;

let == = = ;

% added

picture oddly ; % evenly already defined

evenly := dashpattern(on  3 off 3) ;
oddly  := dashpattern(off 3 on  3) ;

% not perfect, but useful since it removes redundant points.

vardef mfun_straightened(expr sign, p) =
    save _p_, _q_ ; path _p_, _q_ ;
    _p_ := p ;
    forever :
        _q_ := mfun_do_straightened(sign, _p_) ;
        exitif length(_p_) = length(_q_) ;
        _p_ := _q_ ;
    endfor ;
    _q_
enddef ;

vardef mfun_do_straightened(expr sign, p) =
    if length(p)>2 : % was 1, but straight lines are ok
        save pp ; path pp ;
        pp := point 0 of p ;
        for i=1 upto length(p)-1 :
            if round(point i of p) <> round(point length(pp) of pp) :
                pp := pp -- point i of p ;
            fi ;
        endfor ;
        save n, ok ; numeric n ; boolean ok ;
        n := length(pp) ; ok := false ;
        if n>2 :
            for i=0 upto n : % evt hier ook round
                if unitvector(round(point i of pp - point if i=0 : n else : i-1 fi of pp)) <>
                        sign * unitvector(round(point if i=n : 0 else : i+1 fi of pp - point i of pp)) :
                    if ok :
                        --
                    else :
                        ok := true ;
                    fi point i of pp
                fi
            endfor
            if ok and (cycle p) :
                -- cycle
            fi
        else :
            pp
        fi
    else :
        p
    fi
enddef ;

vardef simplified expr p = (
    reverse mfun_straightened(+1,mfun_straightened(+1,reverse p))
) enddef ;

vardef unspiked expr p = (
    reverse mfun_straightened(-1,mfun_straightened(-1,reverse p))
) enddef ;

% path p ;
% p := (2cm,1cm) -- (2cm,1cm) -- (2cm,1cm) -- (3cm,1cm) --
%      (4cm,1cm) -- (4cm,2cm) -- (4cm,2.5cm) -- (4cm,3cm) --
%      (3cm,3cm) -- (2cm,3cm) -- (1cm,3cm) -- (-1cm,3cm) --
%      .5[(-1cm,3cm),(1cm,1cm)] -- (1cm,1cm) -- cycle ;
%
% p := unitcircle scaled 4cm ;
%
% drawpath p ; drawpoints p ; drawpointlabels p ;
% p := p shifted (4cm,0) ; p := straightened p ;
% drawpath p ; drawpoints p ; drawpointlabels p ;
% p := p shifted (4cm,0) ; p := straightened p ;
% drawpath p ; drawpoints p ; drawpointlabels p ;

% new

path originpath ; originpath := origin -- cycle ;

vardef unitvector primary z =
    if abs z = abs origin : z else : z/abs z fi
enddef;

% also new

% vardef anchored@#(expr p, z) = % maybe use the textext variant
%     p shifted (z + (labxf@#*lrcorner p + labyf@#*ulcorner p + (1-labxf@#-labyf@#)*llcorner p))
% enddef ;

% epsed(1.2345)

vardef epsed (expr e) =
    e if e>0 : + eps elseif e<0 : - eps fi
enddef ;

% handy

def withgray primary g =
    withcolor (g,g,g)
enddef ;

% for metafun

if unknown darkred     : color darkred     ; darkred     := .625(1,0,0) fi ;
if unknown darkgreen   : color darkgreen   ; darkgreen   := .625(0,1,0) fi ;
if unknown darkblue    : color darkblue    ; darkblue    := .625(0,0,1) fi ;
if unknown darkcyan    : color darkcyan    ; darkcyan    := .625(0,1,1) fi ;
if unknown darkmagenta : color darkmagenta ; darkmagenta := .625(1,0,1) fi ;
if unknown darkyellow  : color darkyellow  ; darkyellow  := .625(1,1,0) fi ;
if unknown darkgray    : color darkgray    ; darkgray    := .625(1,1,1) fi ;
if unknown lightgray   : color lightgray   ; lightgray   := .850(1,1,1) fi ;

% an improved plain mp macro

vardef center primary p =
    if pair p :
        p
    else :
        .5[llcorner p, urcorner p]
    fi
enddef;

% new, yet undocumented

vardef rangepath (expr p, d, a) =
    if length p>0 :
        (d*unitvector(direction 0 of p) rotated a) shifted point 0 of p
        -- p --
        (d*unitvector(direction length(p) of p) rotated a) shifted point length(p) of p
    else :
        p
    fi
enddef ;

% under construction

vardef straightpath (expr a, b, method) =
    if (method<1) or (method>6)  :
        (a--b)
    elseif method = 1 :
        (a --
            if xpart a > xpart b :
                if ypart a > ypart b :
                    (xpart b,ypart a) --
                elseif ypart a < ypart b :
                    (xpart a,ypart b) --
                fi
            elseif xpart a < xpart b :
                if ypart a > ypart b :
                    (xpart a,ypart b) --
                elseif ypart a < ypart b :
                    (xpart b,ypart a) --
                fi
            fi
        b)
    elseif method = 3 :
        (a --
            if xpart a > xpart b :
                (xpart b,ypart a) --
            elseif xpart a < xpart b :
                (xpart a,ypart b) --
            fi
        b)
    elseif method = 5 :
        (a --
            if ypart a > ypart b :
                (xpart b,ypart a) --
            elseif ypart a < ypart b :
                (xpart a,ypart b) --
            fi
        b)
    else :
        (reverse straightpath(b,a,method-1))
    fi
enddef ;

% handy for myself

def addbackground text t =
    begingroup ;
    save p, b ; picture p ; path b ;
    b := boundingbox currentpicture ;
    p := currentpicture ; currentpicture := nullpicture ;
    fill b t ;
    setbounds currentpicture to b ;
    addto currentpicture also p ;
    endgroup ;
enddef ;

% makes a (line) into an infinite one (handy for calculating
% intersection points

vardef infinite expr p =
    (-infinity*unitvector(direction 0 of p)
    shifted point 0 of p
    -- p --
    +infinity*unitvector(direction length(p) of p)
        shifted point length(p) of p)
enddef ;

% obscure macros: create var from string and replace - and :
% (needed for process color id's) .. will go away

string mfun_clean_ascii[] ;

def register_dirty_chars(expr str) =
    for i = 0 upto length(str)-1 :
        mfun_clean_ascii[ASCII substring(i,i+1) of str] := "_" ;
    endfor ;
enddef ;

register_dirty_chars("+-*/:;., ") ;

vardef cleanstring (expr s) =
    save ss ; string ss, si ; ss = "" ; save i ;
    for i=0 upto length(s) :
        si := substring(i,i+1) of s ;
        ss := ss & if known mfun_clean_ascii[ASCII si] : mfun_clean_ascii[ASCII si] else : si fi ;
    endfor ;
    ss
enddef ;

vardef asciistring (expr s) =
    save ss ; string ss, si ; ss = "" ; save i ;
    for i=0 upto length(s) :
        si := substring(i,i+1) of s ;
        if (ASCII si >= ASCII "0") and (ASCII si <= ASCII "9") :
            ss := ss & char(scantokens(si) + ASCII "A") ;
        else :
            ss := ss & si ;
        fi ;
    endfor ;
    ss
enddef ;

vardef setunstringed (expr s, v) =
    scantokens(cleanstring(s)) := v ;
enddef ;

vardef getunstringed (expr s) =
    scantokens(cleanstring(s))
enddef ;

vardef unstringed (expr s) =
    expandafter known scantokens(cleanstring(s))
enddef ;

% for david arnold:

% showgrid(-5,10,1cm,-10,10,1cm);

def showgrid (expr MinX, MaxX, DeltaX, MinY, MaxY, DeltaY) = % will move
    begingroup
    save size ; numeric size ; size := 2pt ;
    for x=MinX upto MaxX :
        for y=MinY upto MaxY :
            draw (x*DeltaX, y*DeltaY) withpen pencircle scaled
                if (x mod 5 = 0) and (y mod 5 = 0) :
                    1.5size withcolor .50white
                else :
                       size withcolor .75white
                fi ;
        endfor ;
    endfor ;
    for x=MinX upto MaxX:
        label.bot(textext("\infofont " & decimal x), (x*DeltaX,-size)) ;
    endfor ;
    for y=MinY upto MaxY:
        label.lft(textext("\infofont " & decimal y), (-size,y*DeltaY)) ;
    endfor ;
    endgroup
enddef;

% new, handy for:
%
% \startuseMPgraphic{map}{n}
%   \includeMPgraphic{map:germany} ;
%   c_phantom (\MPvar{n}<1) (
%     fill map_germany withcolor \MPcolor{lightgray} ;
%     draw map_germany withpen pencircle scaled 1pt withcolor \MPcolor{darkgray} ;
%   ) ;
%   \includeMPgraphic{map:austria} ;
%   c_phantom (\MPvar{n}<2) (
%     fill map_austria withcolor \MPcolor{lightgray} ;
%     draw map_austria withpen pencircle scaled 1pt withcolor \MPcolor{darkgray} ;
%   ) ;
%   c_phantom (\MPvar{n}<3) (
%   \includeMPgraphic{map:swiss} ;
%     fill map_swiss withcolor \MPcolor{lightgray} ;
%     draw map_swiss withpen pencircle scaled 1pt withcolor \MPcolor{darkgray}  ;
%   ) ;
%   c_phantom (\MPvar{n}<4) (
%   \includeMPgraphic{map:luxembourg} ;
%     fill map_luxembourg withcolor \MPcolor{lightgray} ;
%     draw map_luxembourg withpen pencircle scaled 1pt withcolor \MPcolor{darkgray}  ;
%   ) ;
% \stopuseMPgraphic
%
% \useMPgraphic{map}{n=3}

vardef phantom (text t) = % to be checked
    picture _p_ ;
    _p_ := image(t) ;
    addto _p_ also currentpicture ;
    setbounds currentpicture to boundingbox _p_ ;
enddef ;

vardef c_phantom (expr b) (text t) =
    if b :
        picture _p_ ;
        _p_ := image(t) ;
        addto _p_ also currentpicture ;
        setbounds currentpicture to boundingbox _p_ ;
    else :
        t ;
    fi ;
enddef ;

%D Handy:

def break =
    exitif true fi ;
enddef ;

%D New too:

primarydef p xstretched w = (
    p if (bbwidth (p)>0) and (w>0) : xscaled (w/bbwidth (p)) fi
) enddef ;

primarydef p ystretched h = (
    p if (bbheight(p)>0) and (h>0) : yscaled (h/bbheight(p)) fi
) enddef ;

primarydef p snapped s =
    hide (
        if path p :
            forever :
                exitif (bbheight(p) <= s) and (bbwidth(p) <= s) ;
                p := p scaled (1/2) ;
            endfor ;
        elseif numeric p :
            forever :
                exitif p <= s ;
                p := p scaled (1/2) ;
            endfor ;
        fi ;
    )
    p
enddef ;

% vardef somecolor = (1,1,0,0) enddef ;

% fill OverlayBox withcolor (rcomponent somecolor,gcomponent somecolor,bcomponent somecolor) ;
% fill OverlayBox withcolor (ccomponent somecolor,mcomponent somecolor,ycomponent somecolor,bcomponent somecolor) ;

% This could be standard mplib 2 behaviour:

vardef rcomponent expr p = if rgbcolor  p : redpart     p elseif cmykcolor p : 1 - cyanpart    p else : p fi enddef ;
vardef gcomponent expr p = if rgbcolor  p : greenpart   p elseif cmykcolor p : 1 - magentapart p else : p fi enddef ;
vardef bcomponent expr p = if rgbcolor  p : bluepart    p elseif cmykcolor p : 1 - yellowpart  p else : p fi enddef ;
vardef ccomponent expr p = if cmykcolor p : cyanpart    p elseif rgbcolor  p : 1 - redpart     p else : p fi enddef ;
vardef mcomponent expr p = if cmykcolor p : magentapart p elseif rgbcolor  p : 1 - greenpart   p else : p fi enddef ;
vardef ycomponent expr p = if cmykcolor p : yellowpart  p elseif rgbcolor  p : 1 - bluepart    p else : p fi enddef ;
vardef bcomponent expr p = if cmykcolor p : blackpart   p elseif rgbcolor  p :                 0 else : p fi enddef ;

% draw image       (...) ... ; % prescripts prepended to first, postscripts appended to last
% draw decorated   (...) ... ; % prescripts prepended to each,  postscripts appended to each
% draw redecorated (...) ... ; % prescripts assigned  to each,  postscripts assigned to each
% draw undecorated (...) ... ; % following properties are ignored, existing properties are kept
%
% draw decorated (
%     draw fullcircle scaled 20cm withpen pencircle scaled 20mm withcolor red   withtransparency (1,.40) ;
%     draw fullcircle scaled 15cm withpen pencircle scaled 15mm withcolor green withtransparency (1,.30) ;
%     draw fullcircle scaled 10cm withpen pencircle scaled 10mm withcolor blue  withtransparency (1,.20) ;
% )
%     withcolor blue
%     withtransparency (1,.125) % selectively applied
%     withpen pencircle scaled 10mm
% ;

% vardef image (text imagedata) = % already defined
%     save currentpicture ;
%     picture currentpicture ;
%     currentpicture := nullpicture ;
%     imagedata ;
%     currentpicture
% enddef ;

vardef undecorated (text imagedata) text decoration =
    save currentpicture ;
    picture currentpicture ;
    currentpicture := nullpicture ;
    imagedata ;
    currentpicture
enddef ;


% if metapostversion < 1.770  :
%
%     vardef decorated (text imagedata) text decoration =
%         save mfun_decorated_path, currentpicture ;
%         picture mfun_decorated_path, currentpicture ;
%         currentpicture := nullpicture ;
%         imagedata ;
%         mfun_decorated_path := currentpicture ;
%         currentpicture := nullpicture ;
%         for i within mfun_decorated_path :
%             addto currentpicture
%                 if stroked i :
%                     doublepath pathpart i
%                     dashed dashpart i
%                     withpen penpart i
%                     withcolor colorpart i
%                     decoration
%                 elseif filled i :
%                     contour pathpart i
%                     withpen penpart i
%                     withcolor colorpart i
%                     decoration
%                 elseif textual i :
%                     also i
%                     withcolor colorpart i
%                     decoration
%                 else :
%                     also i
%                 fi
%             ;
%         endfor ;
%         currentpicture
%     enddef ;
%
% else:

    vardef decorated (text imagedata) text decoration =
        save mfun_decorated_path, currentpicture ;
        picture mfun_decorated_path, currentpicture ;
        currentpicture := nullpicture ;
        imagedata ;
        mfun_decorated_path := currentpicture ;
        currentpicture := nullpicture ;
        for i within mfun_decorated_path :
            addto currentpicture
                if stroked i :
                    doublepath pathpart i
                    dashed dashpart i
                    withpen penpart i
                    withcolor colorpart i
                    withprescript prescriptpart i
                    withpostscript postscriptpart i
                    decoration
                elseif filled i :
                    contour pathpart i
                    withpen penpart i
                    withcolor colorpart i
                    withprescript prescriptpart i
                    withpostscript postscriptpart i
                    decoration
                elseif textual i :
                    also i
                    withcolor colorpart i
                    withprescript prescriptpart i
                    withpostscript postscriptpart i
                    decoration
                else :
                    also i
                fi
            ;
        endfor ;
        currentpicture
    enddef ;

% fi ;

vardef redecorated (text imagedata) text decoration =
    save mfun_decorated_path, currentpicture ;
    picture mfun_decorated_path, currentpicture ;
    currentpicture := nullpicture ;
    imagedata ;
    mfun_decorated_path := currentpicture ;
    currentpicture := nullpicture ;
    for i within mfun_decorated_path :
        addto currentpicture
            if stroked i :
                doublepath pathpart i
                dashed dashpart i
                withpen penpart i
                decoration
            elseif filled i :
                contour pathpart i
                withpen penpart i
                decoration
            elseif textual i :
                also i
                decoration
            else :
                also i
            fi
        ;
    endfor ;
    currentpicture
enddef ;

% path mfun_bleed_box ;

% primarydef p bleeded d =
%     image (
%         mfun_bleed_box := boundingbox p ;
%         if pair d :
%             draw p xysized (bbwidth(p)+2*xpart d,bbheight(p)+2*ypart d) shifted -d ;
%         else :
%             draw p xysized (bbwidth(p)+2d,bbheight(p)+2d) shifted (-d,-d) ;
%         fi ;
%         setbounds currentpicture to mfun_bleed_box ;
%     )
% enddef ;

%D New helpers:

def beginglyph(expr unicode, width, height, depth) =
    beginfig(unicode) ; % the number is irrelevant
    charcode := unicode ;
    charwd   := width ;
    charht   := height ;
    chardp   := depth ;
enddef ;

def endglyph =
    setbounds currentpicture to (boundingbox unitsquare xscaled charwd yscaled (charht + chardp) shifted (0,-chardp)) ;
    if known charscale :
        currentpicture := currentpicture scaled charscale ;
    fi ;
    endfig ;
enddef ;

%D Dimensions have bever been an issue as traditional MP can't make that large
%D pictures, but with double mode we need a catch:

newinternal maxdimensions ; maxdimensions := 14000 ;

def mfun_apply_max_dimensions = % not a generic helper, we want to protect this one
    if bbwidth currentpicture > maxdimensions :
        currentpicture := currentpicture if bbheight currentpicture  > bbwidth currentpicture : ysized else : xsized fi maxdimensions ;
    elseif bbheight currentpicture  > maxdimensions :
        currentpicture := currentpicture ysized maxdimensions ;
    fi ;
enddef;

extra_endfig := extra_endfig & "mfun_apply_max_dimensions ;" ;

let dump = relax ;

