%
% PURPOSE: Plot tank test data in Matlab.
% AUTHOR:  Rob McEwen
% DATE:  98/8/27
% MODS: 
% 98/11/19:
% 1) Divided V2lbs by 2 to account for Gary's voltage divider.
% 2) The conditioning amp zero is now at 2.5 volts, or 2048 A/D counts.
% Full scale is 0 to 5 volts. Previously, the cond. amp output range 
% was 106.38 lbs/Volt, centered on about -300 mV.

% 
  nameStr = input(' Data File Name: ','s')
%
% Load the data
  if(nameStr) 
    eval(nameStr)
    nameStrLast = nameStr;
  else
    nameStr = nameStrLast;
  end
% Add lables below so that they match the order of the columns in asciiData.
%
% cnts2V = 5/2^11;           %Datel Conversion from AD counts to volts.
  cnts2V = 5/2^12;           %B&B Conversion from AD counts to volts.
  V2lbs  = 1/.00940*2;       %Conversion from Volts to lbs, from cal. data.
  lbs2N  = 4.448;            %Newtons per pound.
%
  c2N    = cnts2V*V2lbs*lbs2N;
  v2N    = V2lbs*lbs2N;

  lbls = { 'time'
           'ADstatus'
           'ADvalue'
           'vel'
           'temp'
           'thrust'
           'cmdThr'
           'torque'
         };
  sz = size(lbls);
  i = 1;
  while( i <= sz(1) ),
    eval([ char(lbls(i))  ' = '  num2str(i) ';']);
    i = i + 1;
  end

%
% First, compute the initial (first 10 seconds) strain gauge bias reading:
  [mm,idx] = min(abs(data(:,time) - 10 ));
  bias = ones(1,idx)*data(1:idx,ADvalue)/idx*c2N;
  sgauge = data(:,ADvalue)*c2N - bias;
  disp([' The bias is ', num2str(bias), ' Newtons, rms.'])
%
% Now, compute the RMS noise in the first 10 seconds.
  Noise = std(sgauge(1:idx));
  disp(' ')
  disp([' The noise is ', num2str(Noise), ' Newtons, rms.'])
%
% Data strings in the plot:  
  textStr = strvcat(myDate,...
            ['The bias is ',num2str(bias),' N'],...
            ['The 0-10 second noise is ',num2str(Noise), 'N, rms'],...
            ['File Name: ',nameStr]);


  figure(1)
  clf
  sz = size(data);
  N = sz(1);
  plot(data(1:N,time),[data(1:N,vel)])
  ylabel('Shaft Speed, RPM')
  xlabel([' Seconds.  10 Sec Noise = ', num2str(Noise), ' Newtons, rms.'])
  grid on
  title([myDate,' File: ',nameStr,' Bias: ', num2str(bias)])
% val = axis;
% text(aval(2)*.025,aval(4)*.90,textStr)
% legend('Commanded','Measured')

  figure(2)
  clf
  sz = size(data);
  N = sz(1);
  plot(data(1:N,time),[data(1:N,thrust) sgauge(1:N)])
  ylabel('Thrust, Newtons')
  xlabel([' Seconds.  10 Sec Noise = ', num2str(Noise), ' Newtons, rms.'])
  grid on
  title([myDate,' File: ',nameStr,' Bias: ', num2str(bias)])
% val = axis;
% text(aval(2)*.025,aval(4)*.90,textStr)
% legend('Commanded','Measured')

  figure(3)
  clf
  sz = size(data);
  N = sz(1);
  plot(data(1:N,vel).*abs(data(1:N,vel)), sgauge)
  ylabel('Torque, Newton-meters')
  xlabel([' (RPM)^2 '])
  grid on
  title([myDate,' File: ',nameStr,' Bias: ', num2str(bias)])
% val = axis;
% text(aval(2)*.025,aval(4)*.90,textStr)
% legend('Commanded','Measured')

  if sz(2) > 7
     figure(4)
     clf
     sz = size(data);
     N = sz(1);
     plot(data(1:N,time),[data(1:N,torque)/100])
     ylabel('Torque, Newton-meters')
     xlabel([' Seconds.  10 Sec Noise = ', num2str(Noise), ' Newtons, rms.'])
     grid on
     title([myDate,' File: ',nameStr,' Bias: ', num2str(bias)])
  end

  doit = input(' Hardcopy? y => .eps, Y => printer ','s');
  if(doit == 'y' | doit == 'Y')
    eval(['print -f1 -depsc ',nameStr,'_1'])
    eval(['print -f2 -depsc ',nameStr,'_2'])
    eval(['print -f3 -depsc ',nameStr,'_3'])
    if sz(2) > 7, eval(['print -f4 -depsc ',nameStr,'_4']), end
    if( doit == 'Y' )
      eval(['!lpr ',nameStr,'_1.eps'])
      eval(['!lpr ',nameStr,'_2.eps'])
      eval(['!lpr ',nameStr,'_3.eps'])
      if sz(2) > 7, eval(['!lpr ',nameStr,'_4.eps']), end
    end
  end

return

  figure(1)
  clf

  subplot(3,1,1), 
  plot(data(:,time),data(:,thrust))
  ylabel('Newtons')
  title(['Cmded Thrust.',' File: ',nameStr,' Bias: ', num2str(bias),' ',myDate])

  subplot(3,1,2), 
  plot(data(:,time),data(:,vel))
  ylabel('RPM')
  title('Propeller Speed')

  subplot(3,1,3), 
  plot(data(:,time),data(:,temp)),
  xlabel('Time, seconds')
  ylabel('Deg C')
  title('Motor Temperature')

%
% Filter that baby:
  [BB,AA] = butter(4,1/3);  %15/3 = 5 Hz cutoff freq; 4th order
  sgf = filter(BB,AA,sgauge);

  figure(3)
  clf
  sz = size(data);
  N = sz(1);
  plot(data(1:N,time),[data(1:N,thrust) sgf(1:N)])
  ylabel('Thrust, Newtons; 5Hz 4pole Butterworth')
  xlabel([' Seconds'])
  grid on
  legend('Commanded','Measured')
  title([myDate,' File: ',nameStr,' Bias: ', num2str(bias)])
% text(aval(2)*.025,aval(4)*.90,strvcat(textStr,'5 Hz 4th Order Butterworth'))

  figure(3)
  clf

  subplot(2,1,1), 
  plot(data(:,time),sgauge );
  ylabel('Newtons')
  title('Strain Gauge')

  subplot(2,1,2), 
  plot(data(:,time),data(:,ADstatus))
  ylabel('Counts')
  title('Datel A/D Board Status')
  xlabel('Time, seconds')
