function [U,V,outT,outTg,tgrid] = get_UV_from_TRANSP(shot,init,model,params)
%% load TRANSP data
load(sprintf('../../../AUG/AUG_TRANSP_data/TRANSP_%d.mat',shot));

% start at first TRANSP point for simulation
tgrid = params.tgrid(params.tgrid>=data.PCUR.t(1)&params.tgrid<=data.PCUR.t(end));

% get output data from TRANSP to compare
outT = TRANSPdata2RAPTORout(model.rgrid.rho,tgrid,data);
outTg = TRANSPdata2RAPTORout(model.rgrid.rhogauss,tgrid,data);

% get input data from TRANSP

PNB = [data.PINJ01.data;data.PINJ02.data;data.PINJ03.data;data.PINJ04.data;
       data.PINJ05.data;data.PINJ06.data;data.PINJ07.data;data.PINJ08.data];

PEC = [data.PECIN1.data;data.PECIN2.data;data.PECIN3.data;data.PECIN4.data;
       data.PECIN5.data;data.PECIN6.data;data.PECIN7.data;data.PECIN8.data];
   
%%
fnyq = 0.5*mean(1./diff(data.PECIN1.t));
[b,a] = butter(1,5/fnyq,'low');

PEC_filt=filter(b,a,PEC')';
%%
pnb_TRA = interp1(data.PINJ01.t,PNB',tgrid)';
pec_TRA = interp1(data.PECIN1.t,PEC_filt',tgrid)';
Ip_TRA = interp1(data.PCUR.t,data.PCUR.data,tgrid);

%% Define input parameters
% Actuators
U = zeros(model.dims.nu,length(tgrid));
% Plasma current
U(1,:) = Ip_TRA;
% EC power
U(model.echcd.uind(1:8),:) = pec_TRA;

% EC deposition location (variable) - determinined from in specting TRANSP results
switch shot
    case 32148
        params.echcd.rdep(5:8) = [0.2;0.1;0.0;0.5];
        params.echcd.wdep(5:8) = [0.15,0.15,0.15,0.15];
        params.echcd.cd_eff(5:8) = [20,0,0,1];
    case 33864
        params.echcd.rdep(5:8) = [0.2;0.2;0.15;0.4];
        params.echcd.wdep(5:8) = [0.15,0.15,0.15,0.15];
        params.echcd.cd_eff(5:8) = [20,0,-10,20];
        
        params.nbhcd.rdep([5,6,7,8]) =  [0.2,0.6,0.6,0.6];
    otherwise
        error('no ec locations defined for this shot')
end
% NBH power
U(model.nbhcd.uind,:) = pnb_TRA;

% remove NaNs and negatives
U(isnan(U)|U<0) = 0;

% Construct kinetic profiles according to exp data
vb = build_kinetic_profiles(model,init);
v = repmat(vb,1,numel(tgrid));

%% Experimental profiles
% Hmode timing
dLH = 0.1; % duration of LH/HL transition
t_LH1 = 1.2;
t_HL1 = 9.0;
% anonymous function for ramps
rampfun = @(t,tmin,ymin,tmax,ymax) max(ymin,min((ymax-ymin)/(tmax-tmin)*(t-tmin),ymax-ymin)+ymin);

v(model.hmode.vind.activation,:) = rampfun(tgrid,t_LH1,0,t_LH1+dLH,1) - ...
                       rampfun(tgrid,t_HL1,0,t_HL1+dLH,1);

% take from TRANSP directly
v(model.ne.vind,:)    =  model.ne.Lamgauss\outTg.ne;
v(model.ni.vind,:)    =  model.ni.Lamgauss\outTg.ni;

V=v;
return