%% 08/02/2012
% Prove di caratterizzazione harvester V25W come da specifica:
% trave montata come cantilever sull'eccitatore
% tensione in uscita dall'harvester misurata al variare della resistenza
% prove a +-1g e +-0.5g
% R=500,1000,1500,2000,3000,4000 ohm
% m=0,2.4,7.8,15.6 grammi

close all
clear all
clc


%% caricamento dati

% STRUTTURA DEI FILE DI INPUT:

% massa [grammi]
% ampiezza [g]
% freq propria [Hz]
% spostamento tip [V] (peak to peak)
% Resistenza [ohm]  -  dV uscita [V] (peak to peak)

% lista nomi file di input
% ed equazioni delle curve datasheet
lista=struct('file',{'0_1g.dat' '0_05g.dat' '2.4_1g.dat' '2.4_05g.dat'... 
             '7.8_1g.dat' '7.8_05g.dat' '15.6_1g.dat' '15.6_05g.dat'},... 
             'funzione',{'-0.0480*x.^2+0.5338*x' '-0.0446*x.^2+0.2897*x'...
                         '-0.0353*x.^2+0.5861*x' '-0.0353*x.^2+0.2687*x'...
                         '-0.0412*x.^2+1.0892*x' '-0.0334*x.^2+0.5256*x'...
                         '-0.0413*x.^2+1.2357*x' '-0.0408*x.^2+0.6919*x'});

x=[0:.05:30];   
figure('position',[100,100,800,600])
% vettore colori
col=['brbrbrbrbrbrbrbr'];
    
%% plot curve datasheet e risultati sperimentali

for i=1:length(lista)
      
    dati=importdata(['dati/' lista(i).file],'\t');
    
    m=dati(1,1);
    amp=dati(2,1);
    f0=dati(3,1);
    spost=dati(4,1);
    
    R=dati(5:end,1);
    dV=dati(5:end,2)/2./sqrt(2); % /2 ampiezza!
    RMS=dV./sqrt(2);
    
    %P=RMS.^2/R*1000;
    P=dV.^2./R*1000;

    %funzione curva datasheet
    f=inline(lista(i).funzione);
    
    subplot(2,2,(i+mod(i,2))/2)
    hold on
    plot(x,f(x),col(i),'linewidth',2)
    plot(dV,P,[col(i) 'o'])
    ylim([0,inf])
    grid on
    titolo=sprintf('m= %.1f g, f_0= %.1f Hz',m,f0);
    title(titolo);
    xlabel('V [V]')
    ylabel('P [mW]')

end

% salvataggio eps
print -painters -depsc grafici_potenze1.eps
print -painters -dpng grafici_potenze1.png