Merge pull request #117 from MingzeDou/mingze_pull

updated reading digital and optitrack
petersenpeter · Oct 23, 2024 · ca4ae7e · ca4ae7e
2 parents b996c4e + cd7e523
commit ca4ae7e
Show file tree

Hide file tree

Showing 2 changed files with 87 additions and 23 deletions.
diff --git a/calc_CellMetrics/loadOpenEphysDigital.m b/calc_CellMetrics/loadOpenEphysDigital.m
@@ -11,48 +11,102 @@
 % By Peter Petersen
 % petersen.peter@gmail.com
 
-% Load TTL data
-% Each path must contain:
-% 1. timestamps.npy
-% 2: channel_states.npy
-% 3: channels.npy
-% 4: full_words.npy
 
-% TTL_paths = {'TTL_2','TTL_4'};
-% TTL_offset = [0,1];
+% Load TTL data
+% Each path must contain either:
+% Legacy format:
+%   1. timestamps.npy
+%   2. channel_states.npy
+%   3. channels.npy
+%   4. full_words.npy
+% New format:
+%   1. timestamps.npy
+%   2. states.npy
+%   3. sample_numbers.npy
+%   4. full_words.npy
 
 TTL_paths = {};
 epochs_startTime = [];
 ephys_t0 = [];
 
+% Determine file format and set paths
 for i = 1:numel(session.epochs)
-    TTL_paths{i} = fullfile(session.epochs{i}.name,'events','Neuropix-PXI-100.0','TTL_1');
+    % Try new format path first
+    newFormatPath = fullfile(session.epochs{i}.name,'events','Neuropix-PXI-100.ProbeA','TTL');
+    legacyFormatPath = fullfile(session.epochs{i}.name,'events','Neuropix-PXI-100.0','TTL_1');
+
+    % Check which path exists and set accordingly
+    if exist(fullfile(session.general.basePath, newFormatPath), 'dir')
+        TTL_paths{i} = newFormatPath;
+        timestampPath = fullfile(session.general.basePath,session.epochs{i}.name,'continuous','Neuropix-PXI-100.ProbeA','timestamps.npy');
+        isNewFormat = true;  %format flag
+    else
+        TTL_paths{i} = legacyFormatPath;
+        timestampPath = fullfile(session.general.basePath,session.epochs{i}.name,'continuous','Neuropix-PXI-100.1','timestamps.npy');
+        isNewFormat = false;  
+    end
+
     epochs_startTime(i) = session.epochs{i}.startTime;
-    temp = readNPY(fullfile(session.general.basePath,session.epochs{i}.name,'continuous','Neuropix-PXI-100.1','timestamps.npy'));
-    ephys_t0(i) = double(temp(1))/session.extracellular.sr;
+    temp = readNPY(timestampPath);
+    if isNewFormat
+        ephys_t0(i) = double(temp(1));  % timestamps in second according to new format of open ephys
+    else
+        ephys_t0(i) = double(temp(1))/session.extracellular.sr;  
+    end
 end
 
+% Initialize output structure
 openephysDig = {};
-openephysDig.timestamps = epochs_startTime(1) + double(readNPY(fullfile(session.general.basePath, TTL_paths{1},'timestamps.npy')))/session.extracellular.sr - ephys_t0(1);
-openephysDig.channel_states = readNPY(fullfile(session.general.basePath,TTL_paths{1},'channel_states.npy'));
-openephysDig.channels = readNPY(fullfile(session.general.basePath, TTL_paths{1},'channels.npy'));
-openephysDig.full_words = readNPY(fullfile(session.general.basePath,TTL_paths{1},'full_words.npy'));
+
+% Load and process first epoch data
+basePath = fullfile(session.general.basePath, TTL_paths{1});
+timestamps = readNPY(fullfile(basePath,'timestamps.npy'));
+if isNewFormat
+    openephysDig.timestamps = epochs_startTime(1) + double(timestamps) - ephys_t0(1);  % timestamps in second
+else
+    openephysDig.timestamps = epochs_startTime(1) + double(timestamps)/session.extracellular.sr - ephys_t0(1);  % Legacy format
+end
+
+% Check which format to use for first epoch
+if exist(fullfile(basePath,'states.npy'), 'file')
+    openephysDig.channel_states = readNPY(fullfile(basePath,'states.npy'));
+    openephysDig.channels = readNPY(fullfile(basePath,'sample_numbers.npy'));
+    openephysDig.full_words = readNPY(fullfile(basePath,'full_words.npy'));
+else
+    openephysDig.channel_states = readNPY(fullfile(basePath,'channel_states.npy'));
+    openephysDig.channels = readNPY(fullfile(basePath,'channels.npy'));
+    openephysDig.full_words = readNPY(fullfile(basePath,'full_words.npy'));
+end
+
 openephysDig.on{1} = double(openephysDig.timestamps(openephysDig.channel_states == 1));
 openephysDig.off{1} = double(openephysDig.timestamps(openephysDig.channel_states == -1));
 
+% Process additional epochs if present
 if length(TTL_paths) > 1
     openephysDig.nTimestampsPrFile(1) = numel(openephysDig.timestamps);
     openephysDig.nOnPrFile(1) = numel(openephysDig.on{1});
     openephysDig.nOffPrFile(1) = numel(openephysDig.off{1});
+
     for i = 2:length(TTL_paths)
-        timestamps = epochs_startTime(i) + double(readNPY(fullfile(session.general.basePath, TTL_paths{i},'timestamps.npy')))/session.extracellular.sr - ephys_t0(i);
-        openephysDig.timestamps = [openephysDig.timestamps; timestamps];
+        basePath = fullfile(session.general.basePath, TTL_paths{i});
+        timestamps = readNPY(fullfile(basePath,'timestamps.npy'));
 
-        channel_states = readNPY(fullfile(session.general.basePath,TTL_paths{i},'channel_states.npy'));
-        openephysDig.channel_states = [openephysDig.channel_states; channel_states];
+        if isNewFormat
+            timestamps = epochs_startTime(i) + double(timestamps) - ephys_t0(i);  % timestamps in second
+        else
+            timestamps = epochs_startTime(i) + double(timestamps)/session.extracellular.sr - ephys_t0(i);  % Legacy format
+        end
+        openephysDig.timestamps = [openephysDig.timestamps; timestamps];
 
-        openephysDig.channels = [openephysDig.channels;readNPY(fullfile(session.general.basePath, TTL_paths{1},'channels.npy'))];
-        openephysDig.full_words = [openephysDig.full_words; readNPY(fullfile(session.general.basePath,TTL_paths{1},'full_words.npy'))];
+        % Load states and channels based on format
+        if exist(fullfile(basePath,'states.npy'), 'file')
+            openephysDig.channel_states = [openephysDig.channel_states; readNPY(fullfile(basePath,'states.npy'))];
+            openephysDig.channels = [openephysDig.channels; readNPY(fullfile(basePath,'sample_numbers.npy'))];
+        else
+            openephysDig.channel_states = [openephysDig.channel_states; readNPY(fullfile(basePath,'channel_states.npy'))];
+            openephysDig.channels = [openephysDig.channels; readNPY(fullfile(basePath,'channels.npy'))];
+        end
+        openephysDig.full_words = [openephysDig.full_words; readNPY(fullfile(basePath,'full_words.npy'))];
         openephysDig.on{1} = [openephysDig.on{1}; double(timestamps(channel_states == 1))];
         openephysDig.off{1} = [openephysDig.off{1}; double(timestamps(channel_states == -1))];
         openephysDig.nTimestampsPrFile(i) = numel(timestamps);
@@ -78,3 +132,4 @@
 
 % Saving data
 saveStruct(openephysDig,'digitalseries','session',session);
+end
diff --git a/calc_CellMetrics/loadOptitrack.m b/calc_CellMetrics/loadOptitrack.m
@@ -94,7 +94,15 @@
 optitrack_temp.RotationType = dataArray{16}(1);
 optitrack_temp.LenghtUnit = dataArray{18}(1);
 optitrack_temp.CoorinateSpace = dataArray{20}(1);
-optitrack_temp.FrameRate = str2double(dataArray{6}{1});
+
+frameRate_old = str2double(dataArray{6}{1});
+if isempty(frameRate_old) || isnan(frameRate_old)
+    % If old format is invalid, use new format (index 8)
+    optitrack_temp.FrameRate = str2double(dataArray{8}{1});
+else
+    % Use old format value if valid
+    optitrack_temp.FrameRate = frameRate_old;
+end
 
 clear dataArray
 clearvars filename formatSpec fileID dataArray header_length;
@@ -166,7 +174,8 @@
 if parameters.plotFig
     fig1 = figure;
     subplot(1,2,1)
-    plot3(position3D(:,1),position3D(:,2),position3D(:,3)), title('Position'), xlabel('X (cm)'), ylabel('Y (cm)'), zlabel('Z (cm)'),axis tight,view(2), hold on
+    % plot3(position3D(:,1),position3D(:,2),position3D(:,3)), title('Position'), xlabel('X (cm)'), ylabel('Y (cm)'), zlabel('Z (cm)'),axis tight, view(2), hold on
+    plot3(position3D(:,1),position3D(:,2),position3D(:,3)), title('Position'), xlabel('X (cm)'), ylabel('Y (cm)'), zlabel('Z (cm)'),axis tight, hold on
     subplot(1,2,2)
     plot3(position3D(:,1),position3D(:,2),animal_speed), hold on
     xlabel('X (cm)'), ylabel('Y (cm)'),zlabel('Speed (cm/s)'), axis tight