function out=JSP2RAPTORout(rho,time,data)
% out = JSP2RAPTORout(rho,time,data,param);
% leave 'rho' or 'time' empty to leave data on original grids

% table of correspondence between RAPTOR out terms and JSP data
% structure

map = {{'q'   ,    'data.Q'};
     {'iota',      '1./data.Q'};
     {'shear'   ,  'data.SH'};
     {'q0'   ,    'data.Q(1,:)'};
     {'qmin'   ,  'min(data.Q, [], 1)'};
     {'psi' ,     'data.PSI1 * 1e-8 * 2 * pi'};
     {'upl' ,     'data.VOLT'};
     {'te' ,       'data.TE'};
     {'ti' ,       'data.TI'};
     {'ne' ,       'data.NE1 * 1e6'};
     {'ni' ,       'data.NI * 1e6'};
%     {'p'  ,      'data.prof.ptot'},...
     {'ze' ,       'data.ZEFF'};
     {'vt' ,       'data.ANGF'};
     {'Ip' ,       'data.CUR'};
%     {'li3' ,     'data.equi.li'},...
%     {'betaN' ,   'data.equi.betan'},...
%     {'beta' ,    '100*data.gene.beta'},...
%     {'jbs',      'data.neo.jboot'},...
     {'signeo',    'data.SIGM'};
     {'chie',      'data.XE/1e4'};
     {'chii',      'data.XI/1e4'};
     {'pbrem',     'data.QBRE * 1e9'};
     {'prad',      '-data.QRAD * 1e9'};
     {'pei',       'data.QTHX * 1e9'};
     {'poh',       'data.QOH * 1e9'};
     {'pnbe',      'data.QNBE * 1e9'};
     {'pnbi',      'data.QNBI * 1e9'};
     {'pece'  ,    'data.QECE * 1e9'};
     {'pec'  ,     'data.QECE * 1e9'};
     {'pici',      'data.QRFI * 1e9'};
     {'pice',      'data.QRFE * 1e9'};
     {'palphai',      'data.QALI * 1e9'};
     {'palphae',      'data.QALE * 1e9'};
     {'jnb',       'data.JZNB * 3.333333333e-6'};
     {'pfastnb',   'data.WNBD / 10 * (2/3)'};
     {'jic',       'data.JZRF * 3.333333333e-6'};
     {'pfastic',   'data.WRFD / 10 * (2/3)'};
     {'jec',       'data.JZEC * 3.333333333e-6'};
     {'jbs',       'data.JZBS * 3.333333333e-6'};
     {'p',         'data.PR * 1e6'};
%     {'jec',      'data.source.fce.j'},...
%     {'paux',     'data.source.totale.el'},...
%     {'jni',      'data.source.totale.j'},...
%     {'jaux',     'data.source.totale.j-data.source.jboot'},...
%     {'joh',      'data.prof.jmoy-data.source.totale.j'},...
     {'jpar',      'data.JP / 3e5'};
     {'palphae',   'data.QALE * 1e9'};
     {'palphai',   'data.QALI * 1e9'};
%     {'palphatot','data.source.fus.el+data.source.fus.ion'},...
%     {'jbs',      'data.source.jboot'},...
     {'sne' ,      'data.SBD1 * 1e6'};
%     {'Ioh',      'data.gene.ipohm'},...
%     {'Ibs',      'data.gene.iboot'},...
%     {'Iaux',     'data.gene.icd'},...
%     {'Ini',      'data.gene.ini'},...
%     {'Iec',      'data.gene.ifce'},...
%     {'Inb',      'data.gene.iidn'},...
%     {'fbs',      'data.gene.iboot./data.gene.ip'},...
%     {'foh',      'data.gene.ipohm./data.gene.ip'},...
%     {'fni',      'data.gene.ini./data.gene.ip'},...
%     {'tauE',     'data.gene.taue'},...
     {'We'  ,      'sum(data.WE * 1e3)'};
     {'Wi'  ,      'sum(data.WI * 1e3)'};
     {'Wth' ,      'sum((data.WE + data.WI)) * 1e3'};
%     {'Wpol',     'data.gene.wbp'},...
%     {'g1'  ,     'data.equi.grho2.*data.equi.vpr.^2'},... <|gradV|^2>
%     {'g2'  ,     'data.equi.grho2r2.*data.equi.vpr.^2'},... <|gradV|^2/R^2>
%     {'g3'  ,     'data.equi.r2i'},... <1/R^2>
%     {'Vp'  ,     'bsxfun(@times,data.equi.rhomax,data.equi.vpr)'},...
     {'Rgeom',     'data.R(1, :) ./ 100'};
     {'rhopol',     'data.XPSQ'};
    };

%% grids - use input grids or keep old ones
t_jsp   = data.TIME;
r1_jsp = data.XVEC1;
r2_jsp = data.XVEC2;

if isempty(rho)
    interprho = false;
    out.rho = r1_jsp;
    out.rho2 = r2_jsp;
else
    interprho = true;
    out.rho = rho;
end

if isempty(time)
    interptime = false;
    out.time = t_jsp;
else
    interptime = true;
    out.time = time;
end

out.it = 1:numel(out.time);
out.ntime = numel(out.time);
out.nrho = numel(out.rho);

% first eval all strings to get good data
%data = evalstrings(data);

%% interpolate and map to RAPTOR grid
for im=1:numel(map)
    mym = map{im};
    
    try
        jspdata = full(eval(mym{2}));
        
        if size(jspdata, 2) ~= length(t_jsp)
            jspdata = 0.*t_jsp;
        end
        
        % interpolate to new t
        if interptime
            data_on_time = interp1(t_jsp,jspdata',out.time,'linear',NaN)';
        else
            data_on_time = jspdata;
        end

        % if profile data
        if size(jspdata,1)==1 % only time
            out.(mym{1}) = data_on_time';
        else % profile data, interpolate on right rho
            if interprho
                switch size(data_on_time, 1)
                    case numel(r1_jsp)
                        out.(mym{1}) = interp1(r1_jsp,data_on_time,out.rho, 'linear', 'extrap');
                    case numel(r2_jsp)
                        out.(mym{1}) = interp1(r2_jsp,data_on_time,out.rho, 'linear', 'extrap');
                    otherwise
                        error('Could not interpret rho grid')
                end     
            else
                out.(mym{1}) = data_on_time;
            end
        end
    catch err
        warning('could not interpolate %s to %s',mym{2}, mym{1});
        rethrow(err)
    end
end

function dataout=evalstrings(data)

dataout = data;

fields = fieldnames(data);
for ifield = 1:numel(fields)
    myfieldname = fields{ifield};
    myfield = data.(myfieldname);
    if ischar(myfield)
        dataout.(myfieldname) = eval(myfield);
    elseif isa(myfield,'struct')
        dataout.(myfieldname) = evalstrings(myfield);
    end
end
return
