%%
%use this file not the other one
clc 
clear all
close all
%%
%%%files%%%%%%

tilted = 1;
if tilted
    traj_file=readtable('Trajectory_interpolated_tilted_lattice_no_diagAzra.csv');
    spike_file=readtable("encoded_TL_t20p3_b0.8pi_3d_testA.csv");
    neurons_tl = [5,6,7,10,14,18,25,26,33,36,39,41,46];
    neurons = neurons_tl;
    thresh_std = 1.5;
else
    traj_file=readtable('Trajectory_interpolated_aligned_lattice_no_diag.csv');
    spike_file=readtable("encoded_AL_t20p3_b2.51_cossin.csv");
    neurons_al = [2, 4, 6, 8, 11, 16, 18, 24, 26, 30, 33, 34, 41] + 1;
    neurons = neurons_al;
    thresh_std = 1.7;
end
%table_data=csvread(traj_file); %stores traj data in table
table_data = table2array(traj_file);
spike_file=table2array(spike_file);%%loads spike data
plot3(table_data(:,1),table_data(:,2),table_data(:,3))

%%
% % %rotation matrix
% rotm = [[1 0 0]; # for testing only
%     [0 1 0];
%     [0 0 1]];
%x_ang = 0; y_ang = 0;
%tilted = 1;

if tilted
%     x_ang = 45; y_ang = 55;
    rotm = [[ 0.57357644,  0.,          0.81915204],
        [ 0.57922797,  0.70710678, -0.40557979],
        [-0.57922797,  0.70710678,  0.40557979]];
table_data = rotm * table_data';
table_data = table_data';
end
 

%%

for neuron = neurons_tl
%for neuron = 8
%thresh_spike = 0.5;

% for neuron=setdiff(neuron_start:neuron_end,avoid_neurons)
%     neuron
% for neuron = 1:2
    spike_data = spike_file(:,neuron);
      ot_mean = mean(spike_data);
      ot_std = std(spike_data);
      thresh = ot_mean + thresh_std*ot_std;
       firr = spike_data > thresh;
      firposgrid = table_data(firr,:);
      figure
      plot3(table_data(1:50000,1),table_data(1:50000,2),table_data(1:50000,3),'color',[0 0 1 0.1], 'linewidth', 0.001);
      hold on
%       alpha(tt, 0);  
      scatter3(firposgrid(:,1), firposgrid(:,2), firposgrid(:,3),'filled');
%       xlabel('x');
%       ylabel('y');
%       zlabel('z');
%       set(gca, 'YDir','reverse')

%              V1 = [1 1 1];
%              V2 = [6 1 1];
%              V3 = [6 1 6];
%              V4 = [1 1 6];
%              V5 = [1 6 6];
%              V6 = [1 6 1];
%              V7 = [1 6 6];
%              V8 = [6 6 6];
%              V9 = [6 6 1];
%              V10 = [1 6 1];
%              V11 = [1 1 1];
%              V12 = [1 1 6];
%              V13 = [6 1 6];
%              V14 = [6 6 6];
%              V15 = [6 6 1];
%              V16 = [6 1 1];
%              C = [V1;V2;V3;V4;V5;V6;V7;V8;V9;V10;V11;V12;V13;V14;V15;V16];
% %              E = rotm*C';
% %              E = E';
%                E = C;
                
%                  pp = plot3(E(:,1),E(:,2),E(:,3),'black');
% %              pp = gca;
%              rotate(pp, [1 0 0], x_ang);
%              rotate(pp, [0 1 0], y_ang);
% %              xlabel('X');
% %              ylabel('Y');
% %              zlabel('Z');
%                hold on
%              
%             scale=-6:1:6; %scale of reference axes (tilted)
%             %axis A
%             i=1;
%             j=0;
%             k=0;
%             or_sh = repmat(3.5, length(scale),1);
%            
%             i=i*scale;
%             j=j*scale;
%             k=k*scale;                       
%             x_axis = [i; j; k];
%             x_axis = rotm * x_axis;            
%             x_axis(1,:) = x_axis(1,:) + 0.9720;
%             x_axis(2,:) = x_axis(2,:) + 0.8595;
%             x_axis(3,:) = x_axis(3,:) - 5.9217;
%             i = x_axis(1,:);
%             j = x_axis(2,:);
%             k = x_axis(3,:);
%             
%             aa = plot3(i,j,k, 'LineWidth', 1.5,'Color', 'r');
%             legend('A')
%             rotate(aa, [1 0 0], x_ang);
%             rotate(aa, [0 1 0], y_ang);
            hold on
            xlabel('x','Fontsize',15);
            ylabel('y','Fontsize',15);
            zlabel('z','Fontsize',15);
            legend('trajectory', 'firing')
            set(gca, 'FontSize', 15);
            %axis image;
            daspect([1 1 1]);
            %view(30, 30);
            pause();
            savefig(['FF_TL', num2str(neuron)]);
            saveas(gcf, ['FF_TL', num2str(neuron), '.png'])
end
%             
%             
%             
%             %axis B
% %             
%             i = 0;
%             j = 1;
%             k = 0;
%             
%             i=i*scale + or_sh';
%             j=j*scale + or_sh';
%             k=k*scale + or_sh';
%             y_axis = [i; j; k];
% %             y_axis = rotm * y_axis;
%             i = y_axis(1,:);
%             j = y_axis(2,:);
%             k = y_axis(3,:);
%             aa = plot3(i,j,k,'LineWidth', 1.5,'Color', 'g');
%             rotate(aa, [1 0 0], x_ang);
%             rotate(aa, [0 1 0], y_ang);
% %             legend('B')
%             hold on
% %             
%             
%             %axis C
%             i = 0;
%             j = 0;
%             k = 1;
%             
%             i=i*scale + or_sh';
%             j=j*scale + or_sh';
%             k=k*scale + or_sh';
%             z_axis = [i; j; k];
% %             z_axis = rotm * z_axis;
%             i = z_axis(1,:);
%             j = z_axis(2,:);
%             k = z_axis(3,:);
%             aa = plot3(i,j,k,'LineWidth', 1.5, 'Color', 'b');
%             rotate(aa, [1 0 0], x_ang);
%             rotate(aa, [0 1 0], y_ang);
%             legend('A','B','C')
%              hold on
%         
% %              lx2 = E(1:16,1);
% %              ly2 = E(1:16,2);
% %              lz2 = E(1:16,3);
% %              plotOUTLINE;
% %              figure
% %             hax = axes('Parent', fig_); 
% %
% %             pp = plot3(E(:,1),E(:,2),E(:,3),'black');
% %              pp = gca;
% %              rotate(pp, [1 0 0], x_ang);
% %              rotate(pp, [0 1 0], y_ang);
%              
%              hold on
% %       pause();
% %       hold on;
% % end       
        
            %axis off;
            %set(gca, 'YDir','reverse')
%% Firing rate maps
if tilted
    load("Tilted_data_props_b0.8pi_std3d1.3_std2d1.5.mat")
    x_ang = 45;
    y_ang = 55;
else
    load("Aligned_data_props_b2.51_t1.7_2Dt1.5.mat")
    x_ang = 0;
    y_ang = 0;
end


for neuron_firr = neurons_al
    neuron_firr
    % tilted = 0;
        firr_arr=FIRR_3dprops(neuron_firr).firr;
        firr_arr(isnan(firr_arr))=0; 
        col=vals2colormap(firr_arr,'jet',[0.2 1]);        
        if tilted == 0
        figure     
        plotcube_original([20 20 20],[0 0 0],0,[1 1 1])
        pp = gca;
        hold on    
        end
         if tilted == 1
             figure
             V1 = [0 0 0];
             V2 = [20 0 0];
             V3 = [20 0 20];
             V4 = [0 0 20];
             V5 = [0 20 20];
             V6 = [0 20 0];
             V7 = [0 20 20];
             V8 = [20 20 20];
             V9 = [20 20 0];
             V10 = [0 20 0];
             V11 = [0 0 0];
             V12 = [0 0 20];
             V13 = [20 0 20];
             V14 = [20 20 20];
             V15 = [20 20 0];
             V16 = [20 0 0];
             C = [V1;V2;V3;V4;V5;V6;V7;V8;V9;V10;V11;V12;V13;V14;V15;V16];
             E = C;
            pp = plot3(E(:,1),E(:,2),E(:,3),'black');
             rotate(pp, [1 0 0], x_ang);
             rotate(pp, [0 1 0], y_ang);
             hold on           
            scale=-20:1:20; %scale of reference axes (tilted)
            %axis A
%             
            i=1;
            j=0;
            k=0;
            or_sh = repmat(10, length(scale),1);
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            x_axis = [i; j; k];
%             x_axis = rotm * x_axis;
            i = x_axis(1,:);
            j = x_axis(2,:);
            k = x_axis(3,:);     
            aa = plot3(i,j,k, 'LineWidth', 1.5,'Color', 'r');
%             legend('A')
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
            hold on          
            %axis B
%             
            i = 0;
            j = 1;
            k = 0;         
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            y_axis = [i; j; k];
%             y_axis = rotm * y_axis;
            i = y_axis(1,:);
            j = y_axis(2,:);
            k = y_axis(3,:);
            aa = plot3(i,j,k,'LineWidth', 1.5,'Color', 'g');
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
%             legend('B')
            hold on
%                     
            %axis C
            i = 0;
            j = 0;
            k = 1;          
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            z_axis = [i; j; k];
%             z_axis = rotm * z_axis;
            i = z_axis(1,:);
            j = z_axis(2,:);
            k = z_axis(3,:);
            aa = plot3(i,j,k,'LineWidth', 1.5, 'Color', 'b');
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
            %legend('A','B','C')
            hold on
         end     
% %uncomment from here
%             firr_vox = [];
% %             firr_vox = reshape(firr_arr, 20,20,20);
% %             firr_vox = firr_vox';
%             for i = 1:20               
%                 kk = firr_arr(((i-1)*400)+1:(i)*400);
%                 temp = reshape(kk, 20, 20);
%                 firr_vox(:,:,i) = temp;
%             end
% %               gca = rotate(gca, [1 0 0], 45);
% %               gca = rotate(gca, [0 1 0], 55);
% %               gca = rotate(gca, [0 0 1], 0);
%             if tilted == 1
%                 h = vol3d_tilted('Cdata', firr_vox,'texture','3D', 'parent', gca);
%             end
%             if tilted == 0
%                 h = vol3d_aligned('Cdata', firr_vox,'texture','3D', 'parent', gca);
%             end
% %             for i=1:length(h.handles)
% %                 rotate(h.handles(i), [1 0 0], 0, [-0.972, -0.8595, 5.9217]);
% %             end
%             hold on
% %             scatter3(4.26424,24.1421,18.1915,20,'r','filled')
% %             scatter3(4.15121,10,-6.30312,20,'r','filled')
% %             scatter3(-7.32032,10,10.0793,20,'r','filled')
% %             hold on
%             colorbar;
            axis image;
            xlabel('x (voxels)','Fontsize',15);
            ylabel('y (voxels)','Fontsize',15);
            zlabel('z (voxels)','Fontsize',15);
%             set(gca, 'FontSize', 15);
%             axis equal;
%             view(30,30)
% %              rotate(gca, [1 0 0], -90);
% %               rotate(gca, [0 1 0], 225);
% %               rotate(gca, [0 0 1], 225);
% 
% %             set(gca, 'XDir','reverse')
% %             set(gca, 'YDir','reverse')
% %             set(gca, 'ZDir','reverse')
% %             savefig(['FRM_TL', num2str(neuron_firr)]);
% %             saveas(gcf, ['FRM_TL', num2str(neuron_firr), '.png'])
% 
%             pause();
%             %close();
end
% %              axis off;

%% Firing rate maps
if tilted
    load("Tilted_data_props_b0.8pi_std3d1.5_std2d1.5.mat")
    x_ang = 45;
    y_ang = 55;
else
    load("Aligned_data_props_b2.51_t1.7_2Dt1.5.mat")
    x_ang = 0;
    y_ang = 0;
end


for neuron_firr = neurons_tl
    neuron_firr
    % tilted = 0;
        firr_arr=FIRR_3dprops(neuron_firr).firr;
        firr_arr(isnan(firr_arr))=0; 
        col=vals2colormap(firr_arr,'jet',[0.2 1]);        
        if tilted == 0
        figure     
        plotcube_original([20 20 20],[0 0 0],0,[1 1 1])
        pp = gca;
        hold on    
        end
         if tilted == 1
             figure
             V1 = [0 0 0];
             V2 = [20 0 0];
             V3 = [20 0 20];
             V4 = [0 0 20];
             V5 = [0 20 20];
             V6 = [0 20 0];
             V7 = [0 20 20];
             V8 = [20 20 20];
             V9 = [20 20 0];
             V10 = [0 20 0];
             V11 = [0 0 0];
             V12 = [0 0 20];
             V13 = [20 0 20];
             V14 = [20 20 20];
             V15 = [20 20 0];
             V16 = [20 0 0];
             C = [V1;V2;V3;V4;V5;V6;V7;V8;V9;V10;V11;V12;V13;V14;V15;V16];
             E = C;
            pp = plot3(E(:,1),E(:,2),E(:,3),'black');
             rotate(pp, [1 0 0], x_ang);
             rotate(pp, [0 1 0], y_ang);
             hold on           
            scale=-20:1:20; %scale of reference axes (tilted)
            %axis A
%             
            i=1;
            j=0;
            k=0;
            or_sh = repmat(10, length(scale),1);
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            x_axis = [i; j; k];
%             x_axis = rotm * x_axis;
            i = x_axis(1,:);
            j = x_axis(2,:);
            k = x_axis(3,:);     
            aa = plot3(i,j,k, 'LineWidth', 1.5,'Color', 'r');
%             legend('A')
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
            hold on          
            %axis B
%             
            i = 0;
            j = 1;
            k = 0;         
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            y_axis = [i; j; k];
%             y_axis = rotm * y_axis;
            i = y_axis(1,:);
            j = y_axis(2,:);
            k = y_axis(3,:);
            aa = plot3(i,j,k,'LineWidth', 1.5,'Color', 'g');
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
%             legend('B')
            hold on
%                     
            %axis C
            i = 0;
            j = 0;
            k = 1;          
            i=i*scale + or_sh';
            j=j*scale + or_sh';
            k=k*scale + or_sh';
            z_axis = [i; j; k];
%             z_axis = rotm * z_axis;
            i = z_axis(1,:);
            j = z_axis(2,:);
            k = z_axis(3,:);
            aa = plot3(i,j,k,'LineWidth', 1.5, 'Color', 'b');
            rotate(aa, [1 0 0], x_ang);
            rotate(aa, [0 1 0], y_ang);
            %legend('A','B','C')
            hold on
         end     
%uncomment from here
            firr_vox = [];
%             firr_vox = reshape(firr_arr, 20,20,20);
%             firr_vox = firr_vox';
            for i = 1:20               
                kk = firr_arr(((i-1)*400)+1:(i)*400);
                temp = reshape(kk, 20, 20);
                firr_vox(:,:,i) = temp';
            end
            % Rotation matrix around z-axis
%             angle = -90; % Rotate 90 degrees counterclockwise  
%             Rz = [[cosd(angle), -sind(angle), 0],
%                    [sind(angle), cosd(angle), 0],
%                    [0, 0, 1]];
% 
%             % Rotate data             
%             XYZ = firr_vox * Rz.';  
%             Xnew = reshape(XYZ(:,1),size(X));
%             Ynew = reshape(XYZ(:,2),size(Y));
%             Znew = reshape(XYZ(:,3),size(Z));
%             surf(Xnew,Ynew,Znew)
%               gca = rotate(gca, [1 0 0], 45);
%               gca = rotate(gca, [0 1 0], 55);
              rotate(gca, [0 0 1], -90);
            if tilted == 1
                h = vol3d_tilted('Cdata', firr_vox,'texture','3D', 'parent', gca);
%                 h = vol3d_aligned('Cdata', firr_vox,'texture','3D', 'parent', gca);
            end
            if tilted == 0
                h = vol3d_aligned('Cdata', firr_vox,'texture','3D', 'parent', gca);
            end
%             for i=1:length(h.handles)
%                 rotate(h.handles(i), [1 0 0], 0, [-0.972, -0.8595, 5.9217]);
%             end
            hold on
%             scatter3(4.26424,24.1421,18.1915,20,'r','filled')
%             scatter3(4.15121,10,-6.30312,20,'r','filled')
%             scatter3(-7.32032,10,10.0793,20,'r','filled')
%             hold on
            colorbar;
            axis image;
            xlabel('x (voxels)','Fontsize',15);
            ylabel('y (voxels)','Fontsize',15);
            zlabel('z (voxels)','Fontsize',15);
            set(gca, 'FontSize', 15);
            axis equal;
%             view(270, 270);
%              rotate(gca, [1 0 0], -90);
%               rotate(gca, [0 1 0], 180);
%               rotate(gca, [0 0 1], 225);

%             set(gca, 'XDir','reverse')
%             set(gca, 'YDir','reverse')
%             set(gca, 'ZDir','reverse')
            savefig(['FRM_TL', num2str(neuron_firr)]);
%             saveas(gcf, ['FRM_TL', num2str(neuron_firr), '.png'])

%             pause();
            %close();
end
%              axis off;
%%
rotm = [[ 0.70710678, -0.57922797, -0.40557979],
       [ 0.        ,  0.57357644, -0.81915204],
       [ 0.70710678,  0.57922797,  0.40557979]];