% Jon Newman
% Ga Tech

% 16 Apr 2010

% Takes data that has already had PSC's detected and runs a PSC triggered
% average (PSCTA) on the stimulus waveform, if there is one. It does
% separate PSCTA's for EPSC's and IPSC's.

% PSCTA Parameters
% *************************************************************************
dirID = 'C:\TempData\20100409\20100409-1\';
dirdat = dir(strcat([dirID '*.mat']));
causwind = 50; % causal window in ms
acauswind = 10; % acausal window in ms
% *************************************************************************

for j = 1:size(dirdat,1)
    fid = dirdat(j).name;
    load(strcat([dirID fid]));
    
    if isfield(data,'psc')
        if data.stimpara.amp == 0;
            data.psc.peripscsnip = nan;
            data.psc.epscta = nan;
            data.psc.ipscta = nan;
            data.psc.peripsctime = nan;
        else
            % Transform averaging windows into samples
            causwindsamp = floor(data.sampfreq*(causwind/1000));
            acauswindsamp = floor(data.sampfreq*(acauswind/1000));
            
            % Create space for storing peri-epsc stimulus waveform
            data.psc.peripscsnip = cell(1,size(data.psc.coords,2));
            
            for i = 1:size(data.psc.coords,2)
                pscloc = data.psc.coords{i};
                
                % Create peri-psc stimulus waveform data and append to the
                % data vector
                if ~isempty(pscloc)
                    
                    tempsnip = [];
                    epscloc = pscloc(pscloc(:,2)<0);
                    ipscloc = pscloc(pscloc(:,2)>0);
                    
                    tempstim = data.stimulus.vapp(:,i);
                    for k = 1:size(pscloc,1)
                        snipind = [pscloc(k)-causwindsamp pscloc(k)+acauswindsamp];
                        
                        % Get peri-psc stimulus waveform for those psc's
                        % where the peri-psc waveform is not cutoff by the
                        % end of the record
                        if snipind(1) >= 1 && snipind(2) <= length(tempstim)
                            tempsnip = [tempsnip;tempstim(snipind(1):snipind(2))'];
                        else
                            tempsnip = [tempsnip;NaN(1,snipind(2)-snipind(1)+1)];
                        end
                    end
                    data.psc.peripscsnip{i} = tempsnip;
                    data.psc.peripsctime = ((-causwindsamp):acauswindsamp)/(data.sampfreq);
                end
            end
            
            %Find EPSCTA and IPSCTA
            peripscstack = [];
            pscampstack = [];
            for i = 1:size(data.psc.coords,2);
                if ~isempty(data.psc.coords{i})
                    peripscstack = [peripscstack; data.psc.peripscsnip{i}];
                    pscampstack = [pscampstack; data.psc.coords{i}(:,2)];
                end
            end
            
            estack = peripscstack(pscampstack<0,:);
            if ~isempty(estack)
                nans =isnan(estack);
                indnan = find(nans(:,1) ==1);
                estack(indnan,:) = [];
                data.psc.epscta.wave = mean(estack,1);
                if size(estack,1) > 1
                    data.psc.epscta.sd = std(estack);
                else
                    data.psc.epscta.sd = NaN;
                end
            else
                data.psc.epscta = [];
            end
            
            istack = peripscstack(pscampstack>0,:);
            if ~isempty(istack)
                nans =isnan(istack);
                indnan = find(nans(:,1) ==1);
                istack(indnan,:) = [];
                data.psc.ipscta.wave = mean(istack,1);
                if size(istack,1) > 1
                    data.psc.ipscta.sd = std(istack);
                else
                    data.psc.ipscta.sd = NaN;
                end
            else
                data.psc.ipscta = [];
            end
        end
        
        disp(strcat(['Calculated psc-triggered stimulus average. Re-saving file as: ' fid]));
        save(strcat([dirID fid]),'data');
    else
        error('You have not yet found PSCs for these data. Please do this first.');
    end
    
end