`crop/plotarray`:= proc(P,u1,u2,l1,l2) option `Maple Advisor Database 1.01 for Maple 6`, `Copyright (c) 2000 by Robert B. Israel. All rights reserved`; `plot/plotarray`(map(`crop/setview`,P),u1,u2,l1,l2); end; `plots/display_array`:= subs( `plot/plotarray`=`crop/plotarray`, eval(`plots/display_array`)); `crop/setview`:= proc(P) # if P is 2D and has a VIEW option, crop instead local i, PP, r1, r2; if type(P,specfunc(anything,PLOT3D)) then return P fi; for i from 1 to nops(P) do if typematch(op(i,P),VIEW(r1::anything,r2::anything)) then return crop(subsop(i=NULL,P),[`crop/fixdef`(r1), `crop/fixdef`(r2)], true); fi od; P end; `crop/fixdef`:= proc(r) if r = DEFAULT then -infinity..infinity else r fi end; crop := proc(P::specfunc(anything, PLOT), r::[realcons .. realcons, realcons .. realcons]) # P a 2d plot structure # r of the form [a..b,c..d] # optional third argument makes sure there are # points at top left and bottom right if nargs = 2 or hastype(r,infinity) then map(`crop/object`,P,r) else PLOT(POINTS([op([1,1],r),op([2,2],r)],[op([1,2],r),op([2,1],r)], SYMBOL(POINT)), op(map(`crop/object`,P,r))) fi end; `crop/object` := proc(x, r) if type(x, specfunc(anything, CURVES)) then map(`crop/curve`, x, r) elif type(x, specfunc(anything, POLYGONS)) then map(`crop/polygon`, x, r) elif type(x, specfunc(anything, POINTS)) then map(`crop/point`, x, r) elif type(x, specfunc(anything, TEXT)) then if `crop/checkpt`(op(1,x),r) then x else NULL fi else x fi end; `crop/point` := proc(p, r) if type(p, [realcons, realcons]) then if evalb(op([1, 1], r) <= p[1] and p[1] <= op([1, 2], r) and op([2, 1], r) <= p[2] and p[2] <= op([2, 2], r)) then p else NULL fi else p fi end; `crop/checkpt` := proc(p, r) evalb(op([1, 1], r) <= p[1] and p[1] <= op([1, 2], r) and op([2, 1], r) <= p[2] and p[2] <= op([2, 2], r)) end; `crop/curve` := proc(p, r) local res, currentpts, lastpt, i, n, lastin, ptin; if type(p, hfarray) then RETURN(`crop/curve`(convert(p, listlist))) elif not type(p, list([realcons, realcons])) then RETURN(p) fi; n := nops(p); res := NULL; lastpt := p[1]; lastin := `crop/checkpt`(lastpt, r); if lastin then currentpts := lastpt else currentpts := NULL fi; for i from 2 to n do ptin := `crop/checkpt`(p[i], r); if ptin then if lastin then currentpts := currentpts, p[i] else currentpts := `crop/interp`(lastpt, p[i], r), p[i] fi elif lastin then res := res, [currentpts, `crop/interp`(p[i], lastpt, r)] else res := res, `crop/interp`(lastpt, p[i], r) fi; lastpt := p[i]; lastin := ptin od; if currentpts <> NULL then res := res, [currentpts] fi; [res] end; `crop/interp` := proc(p1, p2, r) local v, t1, t2, s1, s2, r1, r2; if p1[1] = p2[1] then if op([1, 1], r) <= p1[1] and p1[1] <= op([1, 2], r) then s1 := 0; s2 := 1 else RETURN(NULL) fi else v := 1/(p2[1] - p1[1]); t1 := (op([1, 1], r) - p1[1])*v; t2 := (op([1, 2], r) - p1[1])*v; if t2 < t1 then s1:= t2; s2:= t1 else s1:= t1; s2:= t2; fi; if s2 < 0 or 1 < s1 then RETURN(NULL) fi; s1 := max(0, s1); s2 := min(1, s2) fi; if p1[2] = p2[2] then if p1[2] < op([2, 1], r) or op([2, 2], r) < p1[2] then RETURN(NULL) fi else v := 1/(p2[2] - p1[2]); t1 := (op([2, 1], r) - p1[2])*v; t2 := (op([2, 2], r) - p1[2])*v; if t2 < t1 then r1:= t2; r2:= t1 else r1:= t1; r2:= t2; fi; s1 := max(s1, r1); s2 := min(s2, r2); if s2 < s1 then RETURN(NULL) fi fi; if s1 = 0 then s2*p2 + (1 - s2)*p1 elif s2 = 1 then s1*p2 + (1 - s1)*p1 else s1*p2 + (1 - s1)*p1, s2*p2 + (1 - s2)*p1 fi end; `crop/polygon` := proc(p, r) local res; if type(p, hfarray) then RETURN(`crop/polygon`(convert(p, listlist))) elif not type(p, list([realcons, realcons])) then RETURN(p) fi; res:= `crop/polygon/halfplane`(p,op([1,1],r),1,true); res:= `crop/polygon/halfplane`(res,op([1,2],r),1,false); res:= `crop/polygon/halfplane`(res,op([2,1],r),2,true); `crop/polygon/halfplane`(res,op([2,2],r),2,false); end; `crop/polygon/halfplane`:= proc(p, a, i, b) # intersect the polygon with x[i] >= a (if b=true) # or x[i] <= a (if b=false), where i=1 or 2 local res, lastpt, lastin, check, current, currin, n, t; if p = [] then RETURN(p) fi; lastpt:= p[1]; if b then check:= t -> evalb(t >= a) else check:= t -> evalb(t <= a) fi; lastin:= check(lastpt[i]); if lastin then res:= lastpt else res:= NULL fi; for n from 2 to nops(p) do current:= p[n]; currin:= check(current[i]); if currin then if lastin then res:= res,current else t:= (current[i]-a)/(current[i]-lastpt[i]); res:= res,(1-t)*current+t*lastpt, current; fi elif lastin then t:= (current[i]-a)/(current[i]-lastpt[i]); res:= res,(1-t)*current+t*lastpt; fi; lastpt:= current; lastin:= currin; od; if lastin <> check(p[1][i]) then t:= (lastpt[i]-a)/(lastpt[i]-p[1][i]); res:= res,(1-t)*lastpt+t*p[1]; fi; [res]; end;