prt_creatertc.m

toolbox of BVQX, This is the access between BV and matlab. It will help you to analysis data from BV

                            rep(rpc).name = sprintf('%s - p1Xp2Xp3', rep(rpc).name);
                    end
                end
            end
            p = rep;
        end
        if numel(p) == 1
            params.params(pc) = p;
        else
            params.params = ...
                [params.params(1:(pc - 1)), p, params.params((pc + 1):end)];
        end
    end
end
numparams = numel(params.params);
parcond = zsz(params.params);
for pc = 1:numparams
    parcond(pc) = params.params(pc).cond;
end

% remove conditions to remove
cond(rcond) = [];
ncond = numel(cond);

% generate one regressor per predictor in MS resolution
nsamp = nvol * prtr;
if ptype
    ncols = 1 + ndnum;
    nrows = nsamp;
else
    ncols = ndcreg;
    nrows = nvol;
end
rtc = zeros(nrows, (ncond + numparams) * ncols);
rtcn = cell(1, (ncond + numparams) * ncols);

% compute first columns per condition / parameter
colcond = ones(1, ncond + 1);
for cc = 2:(ncond + 1)
    colcond(cc) = colcond(cc - 1) + (1 + sum(parcond == (cc - 1))) * ncols;
end

% check size
if (colcond(end) - 1) ~= size(rtc, 2)
    error( ...
        'BVQXfile:InternalError', ...
        'BVQXfile implementation error.' ...
    );
end

% create list of indices for conditions and p

% init RTC
bc2.Type = 'DesignMatrix';
bc2.NrOfDataPoints = nvol;

% create empty cache
ccache = struct;

% iterate over conditions
for cc = 1:ncond
    
    % get first target column
    tcc = colcond(cc);
    
    % any parameters
    pcf = params.params(parcond == cc);

    % get condition on-/offsets and possible parameters
    oo = cond(cc).OnOffsets;
    if ptype
        for dc = [0, nderiv]
            if dc == 0
                rtcn{tcc} = cond(cc).ConditionName{1};
            else
                rtcn{tcc} = ...
                    sprintf('%s - %d.deriv', cond(cc).ConditionName{1}, dc);
            end
            tcc = tcc + 1;
        end
        for pc = 1:numel(pcf)
            pcf(pc).rtccol1 = tcc;
            for dc = [0, nderiv]
                if dc == 0
                    rtcn{tcc} = sprintf('%s x %s', ...
                        cond(cc).ConditionName{1}, pcf(pc).name);
                else
                    rtcn{tcc} = sprintf('%s x %s - %d.deriv', ...
                        cond(cc).ConditionName{1}, pcf(pc).name, dc);
                end
                tcc = tcc + 1;
            end
        end
    else
        namepart = cond(cc).ConditionName{1};
        for nc = 0:(ndcreg - 1)
            rtcn{nc + tcc} = sprintf('%s_D%d', namepart, nc);
        end
        tcc = tcc + ndcreg;
        for pc = 1:numel(pcf)
            for nc = 0:(ndcreg - 1)
                rtcn{nc + tcc} = ...
                    sprintf('%sX%s_D%d', namepart, pcf(pc).name, nc);
            end
            tcc = tcc + ndcreg;
        end
    end

    % iterate over onsets
    for oc = 1:size(oo, 1)
        
        % get first target column again
        tcc = colcond(cc);

        % for volume-PRTs
        if usevol > 0
            oon = (oo(oc, 1) - 1) * usevol + 1 + tshift;
            oof = oo(oc, 2) * usevol + tshift;

        % for millisecond-PRTs
        else
            oon = oo(oc, 1) + 1 + tshift;
            oof = oo(oc, 2) + tshift;
        end
        
        % continue if out of range
        if oon > nsamp || ...
            oof < oon
            continue;
        end
        
        % for HRF
        if ptype
            
            % generate vector of ones
            ov = 1 + oof - oon;
            ccn = sprintf('c%d', ov);

            % already conv
            if isfield(ccache, ccn)
                cv = ccache.(ccn);
            else
                cv = convones(hf, ov);
                cv = cv(:);
                ccache.(ccn) = cv;
            end
            lcv = numel(cv);

            % get indices
            sif = 1;
            sit = lcv;
            if oon < 1
                sif = 2 - oon;
            end
            if (oon + sit - 1) > nsamp
                sit = nsamp + 1 - oon;
            end
            if oon < 1
                sif = sif + 1 - oon;
                oon = 1;
            end

            % add to predictor at onset
            if sit > 0
                
                % iterate over derivatives
                fcc = colcond(cc);
                for dc = [0, nderiv]
                    
                    % no derivative
                    if dc == 0

                        % put HRF into RTC
                        rtc(oon:(oon+sit-sif), tcc) = ...
                            rtc(oon:(oon+sit-sif), tcc) + cv(sif:sit);

                    % derivative
                    else

                        % build derivative and scale to 1000 (match HRF)
                        dcv = [zeros(dc, 1); diff(cv, dc)];
                        dcv = 1000 * dcv / max(abs(dcv));
                        rtc(oon:(oon+sit-sif), tcc) = ...
                            rtc(oon:(oon+sit-sif), tcc) + dcv(sif:sit);
                        
                        % keep track for orthogonalization
                        if ortho
                            olist.(sprintf('r%d', tcc)) = [tcc, (fcc:tcc-1)];
                        end
                    end

                    % increase target column
                    tcc = tcc + 1;
                end
                
                % parameters
                for pc = 1:numel(pcf)
                    
                    % iterate over derivatives
                    for dc = [0, nderiv]

                        % no derivative
                        if dc == 0

                            % put HRF into RTC
                            rtc(oon:(oon+sit-sif), tcc) = ...
                                rtc(oon:(oon+sit-sif), tcc) + ...
                                pcf(pc).pval(oc) * cv(sif:sit);

                            % keep track for orthogonalization
                            if pcf(pc).opts.ortho
                                olist.(sprintf('r%d', tcc)) = [tcc, (fcc:tcc-1)];
                            end
                            
                        % derivative
                        else

                            % build derivative and scale to 1000 (match HRF)
                            dcv = [zeros(dc, 1); diff(cv, dc)];
                            dcv = 1000 * dcv / max(abs(dcv));
                            rtc(oon:(oon+sit-sif), tcc) = ...
                                rtc(oon:(oon+sit-sif), tcc) + ...
                                pcf(pc).pval(oc) * dcv(sif:sit);

                            % keep track for orthogonalization
                            if ortho || pcf(pc).opts.ortho
                                olist.(sprintf('r%d', tcc)) = [tcc, (fcc:tcc-1)];
                            end
                        end

                        % increase target column
                        tcc = tcc + 1;
                    end
                end
            end

        % deconvolution design (possible improvements to be made!)
        else
            for nc = 1:ndcreg
                don = round(oon / prtr) + nc;
                if don > nvol
                    break;
                end
                rtc(don, tcc) = 1;
                tcc = tcc + 1;
            end
            for pc = 1:numel(pcf)
                for nc = 1:ndcreg
                    don = round(oon / prtr) + nc;
                    if don > nvol
                        break;
                    end
                    rtc(don, tcc) = pcf(pc).pval(oc);
                    tcc = tcc + 1;
                end
            end
        end
    end
end

% subsample / normalize RTC (for HRF)
if ptype
    
    % for milliseconds PRTs
    if usevol == 0
        rtc = rtc(1:prtr:end, :);
        
    % for volume-based PRTs
    else
        
        % build new rtc
        nrtc = zeros(nvol, size(rtc, 2));
        
        % iterate over volumes
        for vc = 1:nvol
            nrtc(vc, :) = mean(rtc(1 + (vc - 1) * usevol:vc * usevol, :));
        end
        
        % put back
        rtc = nrtc;
    end

    % normalize predictors
    if params.rnorm ~= 0
        rtc = params.rnorm * rtc / max(rtc(:));
    else
        rtc = rtc / fnorm;
    end
    
    % orthogonalization
    if ortho
        
        % iterate over regressors
        oregs = fieldnames(olist);
        for rc = 1:length(oregs)
            
            % get columns to orthogonalize
            rtcc = olist.(oregs{rc});
            
            % orthoganlize rtcc(1) against the others
            oc = rtcc(1);
            for occ = rtcc(2:end)
                rtc(:, oc) = orthvec(rtc(:, oc), rtc(:, occ));
            end
        end
    end
end

% put into output
bc2.RTCMatrix = rtc;
bc2.PredictorNames = rtcn;
bc2.NrOfPredictors = size(rtc, 2);
bvqxfile_setcont(hfile2.L, bc2);