% File: NMRHourteenPulse.m
%
% Date:   05-Mar-02
% Author: M. Steffen <msteffen@snowmass.stanford.edu
% 
% usage: NMRFourteenPulse(pw90,tau,pulses,templab,phases,expflag)
% 
% pw90 - array of proton and carbon 90 degree pusle widths
%        (eg. pw90=[14.25 9])
% tau - delay time in ms
% pulses - as in the lab manual, array of pulse widths in units of
%          pw90 and phases -0 for x, 1 for y, 2 for mX, 3 for mY]
% templab - which templabl experiment (0,1, or 2)
% phases - values of phase values as determined from the calibration exp.
% expflag - default 0 (do experiment analysis)
%           1 (do everything except for experiment)

function spect = NMRFourPulse(pw90,tau,pulses,templab,phases,expflag) 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Check for correct array size and inputs
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if(nargin~=6)
     expflag=0;
end
if(nargin<4)
     error('NMRFourteenPulses requires 4 inputs!');
end
if(length(pw90)~=2)
     error('pw90 needs to be an array of 2 values!');
end
if( (length(pulses(1,:))~=14) | (length(pulses(:,1))~=2) ) 
     error('pulses needs to be a 2x14 array!');
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CONVERT DELAY TIME TAU INTO SECONDS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

delay=tau/1e3;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SET PULSE SEQUENCE AND DEFINE SOME PARAMETERS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

seqfil='''juniorgrover''';
d1=100;
rcvrs='''yynn''';
sw=400;
mrfb=1000;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SET PARAMETERS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

vjexp(6);
vsetparam('seqfil',seqfil,'%s');

vsetparam('d1',d1,'%d');
vsetparam('rcvrs',rcvrs,'%s');
vsetparam('sw',sw,'%d');
vsetparam('mrfb',mrfb,'%d');

vsetparam('tau',delay,'%d');

vsetparam('ptype',templab,'%d');

vsetparam('pwa',pw90(1),'%d');
vsetparam('pwb',pw90(2),'%d');

vsetparam('pw1a',pulses(1,1),'%d');
vsetparam('ph1a',pulses(1,2),'%d');
vsetparam('pw2a',pulses(1,3),'%d');
vsetparam('ph2a',pulses(1,4),'%d');
vsetparam('pw3a',pulses(1,5),'%d');
vsetparam('ph3a',pulses(1,6),'%d');
vsetparam('pw4a',pulses(1,7),'%d');
vsetparam('ph4a',pulses(1,8),'%d');
vsetparam('pw5a',pulses(1,9),'%d');
vsetparam('ph5a',pulses(1,10),'%d');
vsetparam('pw6a',pulses(1,11),'%d');
vsetparam('ph6a',pulses(1,12),'%d');
vsetparam('pw7a',pulses(1,13),'%d');
vsetparam('ph7a',pulses(1,14),'%d');

vsetparam('pw1b',pulses(2,1),'%d');
vsetparam('ph1b',pulses(2,2),'%d');
vsetparam('pw2b',pulses(2,3),'%d');
vsetparam('ph2b',pulses(2,4),'%d');
vsetparam('pw3b',pulses(2,5),'%d');
vsetparam('ph3b',pulses(2,6),'%d');
vsetparam('pw4b',pulses(2,7),'%d');
vsetparam('ph4b',pulses(2,8),'%d');
vsetparam('pw5b',pulses(2,9),'%d');
vsetparam('ph5b',pulses(2,10),'%d');
vsetparam('pw6b',pulses(2,11),'%d');
vsetparam('ph6b',pulses(2,12),'%d');
vsetparam('pw7b',pulses(2,13),'%d');
vsetparam('ph7b',pulses(2,14),'%d');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% RUN EXPERIMENT AND GET DATA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if(expflag==0)
fprintf(' %%%% starting experiment now. \n');
unix('sleep 1');
vgawait(6);
s=sprintf('vgetdatnoplot(6,%d)',templab);
eval(s);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% READ IN DATA AND FFT USING V2SDAT2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

s=sprintf('sdp%d=v2sdat2j(%d,1,0)',templab,templab); eval(s);
s=sprintf('sdc%d=v2sdat2j(%d,2,0)',templab,templab); eval(s);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CALCULATE PINTS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

s=sprintf('peakfrq_p%d = peaksj(sdp%d,0);',templab,templab); eval(s);
s=sprintf('peakfrq_c%d = peaksj(sdc%d,0);',templab,templab); eval(s);

s=sprintf('peaks_p%d=pintj(sdp%d,peakfrq_p%d,10);',templab,templab,templab); eval(s);
s=sprintf('peaks_c%d=pintj(sdc%d,peakfrq_c%d,10);',templab,templab,templab); eval(s);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ASSIGN SPECT STRUCTURE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

s=sprintf('spect.hfreq =  sdp%d.fdat;',templab); eval(s);
s=sprintf('spect.cfreq =  sdc%d.fdat;',templab); eval(s);
s=sprintf('spect.hspect = sdp%d.spect*exp(-i*phases(1));', templab); eval(s);
s=sprintf('spect.cspect = sdc%d.spect*exp(-i*phases(2));', templab); eval(s);
s=sprintf('spect.hsfo =   sdp%d.sfo;',templab); eval(s);
s=sprintf('spect.csfo =   sdc%d.dfo;',templab); eval(s);
s=sprintf('spect.dt =     sdc%d.dt;',templab); eval(s);
s=sprintf('spect.hpeaks = peaks_p%d*exp(-i*phases(1));',templab); eval(s);
s=sprintf('spect.cpeaks = peaks_c%d*exp(-i*phases(2));',templab); eval(s);
s=sprintf('spect.hphase = phases(1);',templab); eval(s);
s=sprintf('spect.cphase = phases(2);',templab); eval(s);
s=sprintf('spect.pp =     sdc%d.pp;',templab); eval(s);



figure(1)
plot(spect.hfreq,real(spect.hspect));
title('proton spectrum');
ylabel('signal [arb. units]');
xlabel('frequency from hsfo');

figure(2)
plot(spect.cfreq,real(spect.cspect));
title('carbon spectrum');
ylabel('signal [arb. units]');
xlabel('frequency from csfo');
