(* :Author: Christopher Moretti *) 
(* :Summary: This package contains the command 
    RiemannPlot,which creates plots of 
    functions and area estimates. *) 
(* :History: Version 1.0 by Christopher Moretti, 1997 *) 

Clear[RiemannPlot,LeftPoints,RightPoints,MidPoints] 
LeftPoints="l";RightPoints="r";MidPoints="m"; 

Options[RiemannPlot]={HeightPoints->LeftPoints,Intervals->10,Estimate->False, 
    PointPosition->False} 

<<Utilities`FilterOptions` 

RiemannPlot[function_,{x_,a_,b_},opts___]:= 
 Block[{i,dx,n,height,shift,linelist,bareplot,est,lab,pointpos,plotopts}, 
  plotopts= 
      Union[{FilterOptions[Plot,{opts}]},{DisplayFunction->$DisplayFunction}]; 
  
  {height,n,est, 
        pointpos}= {HeightPoints,Intervals,Estimate,PointPosition}/.{ 
            opts}/.Options[RiemannPlot]; 
  If[!(MemberQ[{height},"l"]||MemberQ[{height},"r"]|| MemberQ[{height},"m"]), 
      Return["Invalid HeightPoints choice!"]]; 
  shift=Which[NumberQ[pointpos]&& pointpos>=0  && pointpos<=1,pointpos, 
        height=="l",0,height== "r",1,height=="m",1/2,True,0]; 
  If[NumberQ[pointpos]&&!( pointpos>=0  && pointpos<=1), 
      Print["Invalid point position used - using HeightPoints option \ 
istead"]]; 
  dx=(b-a)/n; 
  linelist=Flatten[Table[ 
   {Line[{ {a+i dx,0},{a+i dx,function/. (x->(a+i dx+shift dx))}}], 
     Line[{ 
                 {a+i dx,function/. (x->(a+i dx+shift dx))},{a+(i+1)dx, 
                  function/. (x->(a+i dx+shift dx))}}], 
     Line[{ {a+(i+1)dx,function/. (x->(a+i dx+shift dx))},{a+(i+1)dx,0} }]}, 
    {i,0,n-1}]]; 
  bareplot=Plot[ function,{x,a,b},DisplayFunction->Identity]; 
  Show[{bareplot,Graphics[linelist]},plotopts]; 
  If[est,Print["The estimated signed area is ", 
    Sum[ dx*(function /.(x->a+i dx +shift dx)),{i,0,n-1}]//N, 
        " square units."]]; 
 ] 
RiemannPlot::usage= 
  "RiemannPlot[f[x],{x,a,b},options] creates a plot of the graph of y=f[x] on \ 
the interval [a,b] and displays rectangles which approximate the signed area \ 
under the curve. The default settings generate 10 rectangles whose heights \ 
are given by the curve at left-hand endpoints.  The number of rectangles can \ 
be changed by the option Intervals (Intervals->10 is default).  The option \ 
HeightPoints (with valid values LeftPoints, RightPoints, Midpoints) changes \ 
where the rectangles get their height.  The option Estimate->True returns a \ 
numerical estimate of the area as well. 
  The option PointPosition->a (a between 0 and 1) overrides Heightpoints to \ 
use the points which are ath of the way across each subinterval. 
  RiemannPlot also accepts all of the options for the Plot command, such as \ 
Ticks, PlotRange, and PlotLabel. 
 "