function area=numtrap(option)
%NUMINT  Shows Trap's rule for a function
%USE:    Type  numtrap  at the MATLAB prompt

%Peter Dunn
%26 April 1999

if nargin==0,

	figure('Name','Trapezoidal Rule', 'tag','tagTrap',...
   		 'NumberTitle','off'	);

   enterframe = uicontrol('Style','Frame','Units','normalized',...
            'Position',[0.02 0.02 0.79, 0.24]);

   funtext = uicontrol('Style','text','Units','normalized',...
             'Position',[0.03 0.12 0.16,0.08],'String','Function:',...
             'FontSize',14,...
             'HorizontalAlignment','center');
   fun = uicontrol('Style','edit','Units','normalized',...
         'Position',[0.19 0.15 0.47 0.08],...
         'FontSize',14,...
         'ToolTip','The function to integrate',...
         'Callback','numtrap(5);',...
         'String','sin(x.^2)','tag','tagTrapFun' );

   lowertext = uicontrol('Style','text','Units','normalized',...
         'FontSize',14,...
         'Position',[0.03 0.03 0.11 0.08],'String','Min x:' );
   lower = uicontrol('Style','edit','Units','normalized',...
         'FontSize',14,...
         'ToolTip','The lower integration limit',...
         'Position',[0.15 0.03 0.08 0.08],...
         'Callback','numtrap(6);',...
         'String','0','tag','tagTrapLower');

   uppertext = uicontrol('Style','text','Units','normalized',...
         'FontSize',14,...
         'Position',[0.24 0.03 0.11 0.08],'String','Max x:');
   upper = uicontrol('Style','edit','Units','normalized',...
         'FontSize',14,...
         'ToolTip','The upper integration limit',...
         'Position',[0.36 0.03 0.08 0.08],...
         'Callback','numtrap(6);',...
         'String','2','tag','tagTrapUpper');

   numptstext = uicontrol('Style','text','Units','normalized',...
         'FontSize',14,...
         'Position',[0.45 0.03 0.12 0.08],'String','Panels:');
   numpts = uicontrol('Style','edit','Units','normalized',...
         'FontSize',14,...
         'ToolTip','Number of panels to use',...
         'Callback','numtrap(10);',...
         'Position',[0.58 0.03 0.08 0.08],'String','6',...
         'tag','tagTrapNum');

   speedtextup = uicontrol('Style','text','Units','normalized',...
         'FontSize',14,...
         'Position',[0.70 0.20 0.10 0.04],'String','50');
   speedtextlo = uicontrol('Style','text','Units','normalized',...
         'FontSize',14,...
         'ToolTip','Changes the number of panels to use',...
         'Position',[0.70 0.03 0.10 0.04],'String','2');
   speed = uicontrol('Style','slider','Units','normalized',...
         'Position',[0.72 0.08 0.06 0.12],'Min',2,'Max',50,...
         'Callback','numtrap(11);',...
         'SliderStep',[1/48 4/48],  'Value',6,'tag','tagTrapSlide');

   gobtn = uicontrol('Style','pushbutton','Units','normalized',...
             'FontSize',14,...
             'ToolTip','Approximates the integral',...
             'Position',[0.83 0.14 0.12 0.12],...
             'String','GO','Callback','numtrap(4); ');
   quitbtn = uicontrol('Style','pushbutton','Units','normalized',...
             'FontSize',14,...
             'ToolTip','Closes the program',...
             'Position',[0.83 0.02 0.12 0.12],...
             'String','QUIT','Callback','delete(gcf);');

   uimenu('Label','DEMONSTRATIONS', 'tag','tagTrapDemo');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(20);',...
           'Label','Difficult Regions: sqrt(x)');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(21);',...
           'Label','Difficult Regions: exp( -abs(x) )');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(22);',...
           'Label','Difficult Regions: 1 ./ ( 1 + x.^6 )');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(23);',...
           'Label','Difficult Regions: besselj(0.2, x)');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(24);',...
           'Label','Difficult Regions: sin( 1 ./ x )');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(25);',...
           'Label','Under Estimation: x.^2 + 4');
   uimenu( findobj('tag', 'tagTrapDemo'), ...
           'Callback','numtrap(26);',...
           'Label','Over Estimation: exp( -x )');

	axes('Position',[0.10 0.35 0.85, 0.6]);
	axis( [ 0 2 0 1 ] );

   set( gca,'FontSize',16);

   zoom;

   %Flag that othing has been printed yet:
   set( findobj('tag','tagTrap'), 'UserData',0);

   %numtrap(4);


elseif nargin==1,

   if option == 3, %NUMBER REGIONS

      np = get( findobj(gcf,'tag','tagTrapNum'),'String');
      np = floor(str2num(np));

      if np<1,
         watchoff;
         herr=errordlg('Need at least two points!','Error!');
         set(herr,'WindowStyle','modal');
         return;
      end;

      if np>50,
         watchoff;
         herr=errordlg(...
               'More than 50 points will slow MATLAB down too much!',...
               'Error!');
         set(herr,'WindowStyle','modal');
         set( findobj(gcf,'tag','tagTrapNum'),'String','50');
         np = 50;
         return;
      end;

   elseif option == 4, %FUNCTION

      prt = get( findobj('tag','tagTrap'),'UserData');
      set( findobj('tag','tagTrap'),'UserData',1);

      if ( prt==0 ),

         %Set up Command Window output
         numtrap(5);

      end;

      watchon;

      delete(gca);

      numtrap(3);

      f = get( findobj(gcf,'tag','tagTrapFun'),'String');
      
      xlow = get( findobj(gcf,'tag','tagTrapLower'), 'String');
      xupp = get( findobj(gcf,'tag','tagTrapUpper'), 'String');
      numpt = get( findobj(gcf,'tag','tagTrapNum'), 'String');

      %Plot the function
      xl = str2num(xlow);
      xu = str2num(xupp);
      np = str2num(numpt);

      axes('Position',[0.10 0.35 0.85, 0.6]);

      x = linspace(xl,xu,1000);

      prob = 0;
      eval(['fy=',f,';'],'fy=NaN;');
      fx = x;

      if any(isnan(fy)),
         watchoff;
         herr=errordlg(...
              ['Couldn''t evaluate the function!  Check the function ',...
              'and the evaluation limits.'],'Error evaluating function');
         set(herr,'WindowStyle','modal');
         return;
      end;

      hold on;

      plot( fx, fy, 'g','LineWidth',4 );
      plot( [min(fx),max(fx)],[0 0],'k-');
      hold on;

      %Determine interval info
      h = ( xu - xl ) / np ;
      xx = [ xl : h : xu ];
      x = xx;

      approxarea = 0;

      for i=0:np-1,
 
         left = xl + h*i;
         right = left+h;
         x = left;
 
         stop = 0;
         ht = eval(f);
 
         x=right;
         htr = eval(f);
 
         [area, hndl] = plottrap( left, right, ht, htr );
         approxarea = approxarea + area;
 
      end;

      %LEGEND is slow and causes the figure to be redrawn... speed up?
      legend( [hndl],['Integral = ',num2str(approxarea)],0);

      set( gca,'FontSize',16);

      fprintf( '%5i   %15.7f\n', np, approxarea);

      watchoff;

   elseif option==5, %Header for output

      %Set up Command Window output

      f = get( findobj('tag','tagTrapFun'),'String');
      lo = get( findobj('tag','tagTrapLower'),'String');
      hi = get( findobj('tag','tagTrapUpper'),'String');
      disp(' ');
      disp( [ 'Integrating ',f,' between ', ...
              num2str(lo),' and ',num2str(hi) ] );
      disp(' ');
      disp('Number of    Approximate');
      disp(' Panels        Integral');
      disp('~~~~~~~~~    ~~~~~~~~~~~~');

   elseif option==6,

      %Flag need for new header

      set( findobj('tag','tagTrap'), 'UserData',0);

   elseif option==10, %Number Panels changed

      val = get( findobj('tag','tagTrapNum'), 'String');
      %Convert to number
      val = val{1};
      val = eval(val);
      val = floor(val);

      if val > 50,
         val = 50;
         set( findobj('tag','tagTrapNum'), 'String','50');
      elseif val <= 0,
         val = 2;
         set( findobj('tag','tagTrapNum'), 'String','2');
      else
         set( findobj('tag','tagTrapNum'), 'String',num2str(val));
      end;

      set( findobj('tag','tagTrapSlide'), 'Value', val );

   elseif option==11, %Slider moves

      val = get( findobj('tag','tagTrapSlide'), 'Value' );
      val = floor(val);
      set( findobj('tag','tagTrapNum'), 'String',num2str(val) );

      numtrap(4);

   elseif option==20,

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','sqrt( x )' );
      set( findobj('tag','tagTrapLower'),...
           'String','0');
      set( findobj('tag','tagTrapUpper'),...
           'String','4');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==21,

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','exp( -abs(x) )' );
      set( findobj('tag','tagTrapLower'),...
           'String','-2');
      set( findobj('tag','tagTrapUpper'),...
           'String','2');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==22,

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','1 ./ ( 1 + x.^6 )' );
      set( findobj('tag','tagTrapLower'),...
           'String','-3');
      set( findobj('tag','tagTrapUpper'),...
           'String','4');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==23,

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','besselj( 0.2, x )' );
      set( findobj('tag','tagTrapLower'),...
           'String','0');
      set( findobj('tag','tagTrapUpper'),...
           'String','1');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==24,

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','sin( 1 ./ x )' );
      set( findobj('tag','tagTrapLower'),...
           'String','0.1');
      set( findobj('tag','tagTrapUpper'),...
           'String','1');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==25, %Under Estimation

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','x.^2 + 4 ' );
      set( findobj('tag','tagTrapLower'),...
           'String','-1');
      set( findobj('tag','tagTrapUpper'),...
           'String','2');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   elseif option==26, %Over Estimation

      numtrap(6);

      set( findobj('tag','tagTrapFun'), ...
           'String','exp( -x )' );
      set( findobj('tag','tagTrapLower'),...
           'String','-4');
      set( findobj('tag','tagTrapUpper'),...
           'String','1');
      set( findobj('tag','tagTrapSlide'),...
           'Value',2);
      set( findobj('tag','tagTrapNum'),...
           'String','2');

      numtrap(4);

   end;

end;

return;


%%%%%%%%%%%%%%%%
%  SUBFUNCTION %
%%%%%%%%%%%%%%%%

 
function [area, h]=plottrap(xlow, xupp, ht, htr)
%PLOTS the trapeziod
 
h = plot( [xlow xlow xupp xupp], [0 ht htr 0] , 'r-', ...
          'LineWidth',2 );
 
area = (xupp - xlow ) * ( htr + ht ) / 2;

plot( xlow, ht, 'r.','MarkerSize',20);
plot( xupp, htr, 'r.','MarkerSize',20);

set(gca, 'FontSize', 14);
 
return;