From 59765563bdffef2467ef6827f7129865c449c1dc Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 10 Apr 2026 14:42:37 +0000 Subject: [PATCH 1/3] Fix reading digital channels in neuropixelsGLX via high-level read() API ndr.reader.read() routed 'digital_in' / 'di' channel types to readevents_epochsamples (intended for event/marker data), but ndr.reader.neuropixelsGLX only implements readchannels_epochsamples for digital_in (the 16-bit NIDQ sync word is regularly sampled at the same rate as the analog channels). The event path falls through to the abstract base implementation, which returns [], so callers like r.read(epoch,'di1','t0',0,'t1',10) got empty data and time. Route digital_in/digital_out through readchannels_epochsamples in the high-level read() switch, the same way analog_in is handled. This preserves the native int16 digital word (it is a 16-bit signal, not 0/1) and mirrors how intan_rhd already treats its digital lines in readchannels_epochsamples. Extend the neuropixelsGLX reader test to exercise the sync channel through r.read('di1',...) with a 16-bit ramp pattern to guard against regressions that would collapse the word to 0/1 or re-break the routing. --- +ndr/+test/+reader/+neuropixelsGLX/test.m | 23 +++++++++++++++++++++++ +ndr/reader.m | 3 ++- 2 files changed, 25 insertions(+), 1 deletion(-) diff --git a/+ndr/+test/+reader/+neuropixelsGLX/test.m b/+ndr/+test/+reader/+neuropixelsGLX/test.m index eb7bdf4..80f7bd3 100644 --- a/+ndr/+test/+reader/+neuropixelsGLX/test.m +++ b/+ndr/+test/+reader/+neuropixelsGLX/test.m @@ -42,6 +42,11 @@ function test(varargin) for c = 1:nNeuralChans data(:, c) = int16(round(500 * sin(2 * pi * c * t_vec))); end + % Sync channel is a 16-bit digital word (packed digital lines), not 0/1. + % Use a ramping pattern that exercises multiple bits to confirm the full + % 16-bit value is preserved end-to-end. + sync_data = int16(mod((0:nSamples-1), 2^15)); + data(:, nTotalChans) = sync_data(:); % Write binary fid = fopen(binfile, 'w', 'ieee-le'); @@ -102,6 +107,24 @@ function test(varargin) disp(['Max error vs expected: ' num2str(max_error)]); assert(max_error == 0, 'Data mismatch detected!'); + % Read the digital (sync) channel through the high-level read() API. + % This is the code path that previously returned [] because + % ndr.reader.read() routed 'digital_in' to readevents_epochsamples + % (which is abstract for this format). + t0 = 0; + t1 = (nSamples-1) / SR; + [d_di, t_di] = r.read({metafile}, 'di1', 't0', t0, 't1', t1); + disp(['Read ' int2str(size(d_di, 1)) ' samples from channel di1 via r.read().']); + assert(~isempty(d_di), 'Digital read returned empty data.'); + assert(~isempty(t_di), 'Digital read returned empty time.'); + assert(size(d_di, 1) == nSamples, ... + sprintf('Digital sample count mismatch: got %d, expected %d.', ... + size(d_di, 1), nSamples)); + % The full 16-bit word must be preserved (not collapsed to 0/1). + di_error = max(abs(double(d_di) - double(sync_data(:)))); + disp(['Max digital error vs expected: ' num2str(di_error)]); + assert(di_error == 0, 'Digital (16-bit) word mismatch detected!'); + disp('All checks passed.'); if plotit diff --git a/+ndr/reader.m b/+ndr/reader.m index e86e694..5d38382 100644 --- a/+ndr/reader.m +++ b/+ndr/reader.m @@ -131,7 +131,8 @@ if b, switch (channelstruct(1).ndr_type), % readchannels_epochsamples - case {'analog_input','analog_output','analog_in','analog_out','ai','ao'}, + case {'analog_input','analog_output','analog_in','analog_out','ai','ao', ... + 'digital_input','digital_output','digital_in','digital_out','di','do'}, if ~useSamples, % must compute the samples to be read s0 = round(1+t0*channelstruct(1).samplerate); s1 = round(1+t1*channelstruct(1).samplerate); From 3b762cd757a4822863bf0924351a74b89598c3c6 Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 10 Apr 2026 14:53:18 +0000 Subject: [PATCH 2/3] Expose neuropixelsGLX digital word as per-bit di channels In NDR each digital_in channel represents a single bit. Previously the SpikeGLX reader exposed the entire packed digital word as a single 'di1' channel, which conflicted with the NDR convention and made it impossible for callers to ask for an individual line. Detect digital width from metadata: - NIDQ: n_digital_lines = 8 * (niXDBytes1 + niXDBytes2). For the user-reported file (niXDBytes1=1) this gives 8 lines (di1..di8). A future file with niXDBytes1=2 would yield 16, and adding port1 via niXDBytes2 extends further. Falls back to 16 * n_dw_chans when the byte fields are absent. - IMEC: n_digital_lines = 16 * n_sync_chans, so the int16 sync word is exposed as di1..di16 (bit 6 is the SMA sync input in practice; other bits are available if a recording uses them). header.m gains two new fields, n_digital_word_cols (number of int16 columns in the file holding digital data, always at the end) and n_digital_lines (number of single-bit lines exposed). Existing fields are unchanged so samples2volts and the NIDQ tests keep working. getchannelsepoch now emits one digital_in entry per bit (di1..di_n_digital_lines) instead of a single di1. readchannels_epochsamples('digital_in', channels, ...) maps each 1-based line index to (DW column, bit position), reads each unique DW column once via read_samples, and extracts the requested bit with bitget. Returns int16 0/1 values, one column per requested line. Out-of-range lines raise ndr:reader:neuropixelsGLX:DigitalLineOutOfRange. Tests: - +ndr/+test/+reader/+neuropixelsGLX/test.m: writes a counter pattern to the IMEC sync column and verifies di1, di8 and di12 return the correct bit via the high-level r.read() API. Adds a NIDQ section that mirrors the user's snsMnMaXaDw=0,0,8,1 + niXDBytes1=1 configuration, asserts exactly 8 di channels are listed, reads di1 via r.read() (the exact failing call from the bug report), and verifies di9 raises DigitalLineOutOfRange. - tools/tests/.../TestNeuropixelsGLX.m: testGetChannelsEpoch updated to expect 16 digital lines after the neural channels instead of a single di1. --- +ndr/+format/+neuropixelsGLX/header.m | 36 ++++- +ndr/+reader/neuropixelsGLX.m | 55 +++++-- +ndr/+test/+reader/+neuropixelsGLX/test.m | 141 +++++++++++++++--- .../+unittest/+reader/TestNeuropixelsGLX.m | 18 ++- 4 files changed, 214 insertions(+), 36 deletions(-) diff --git a/+ndr/+format/+neuropixelsGLX/header.m b/+ndr/+format/+neuropixelsGLX/header.m index 6582a37..0600fde 100644 --- a/+ndr/+format/+neuropixelsGLX/header.m +++ b/+ndr/+format/+neuropixelsGLX/header.m @@ -67,7 +67,13 @@ % Number of saved channels info.n_saved_chans = str2double(meta.nSavedChans); - % Parse snsApLfSy or snsMnMaXaDw to determine neural vs sync channels + % Parse snsApLfSy or snsMnMaXaDw to determine neural vs sync channels. + % Also compute: + % n_digital_word_cols : number of int16 columns in the .bin file + % that hold digital word data (stored last). + % n_digital_lines : number of single-bit digital lines exposed. + % SpikeGLX packs digital data into int16 columns; + % each bit is an independent digital line. if isfield(meta, 'snsApLfSy') % imec stream: AP,LF,SY counts counts = sscanf(meta.snsApLfSy, '%d,%d,%d'); @@ -82,6 +88,12 @@ info.stream_type = 'ap'; info.n_neural_chans = counts(1); end + % IMEC sync word is an int16; each sync column provides 16 bits. + % In practice only bit 6 is the SMA sync input, but all 16 bits + % are exposed as independent digital lines so callers can pick + % whichever they need. + info.n_digital_word_cols = info.n_sync_chans; + info.n_digital_lines = 16 * info.n_sync_chans; elseif isfield(meta, 'snsMnMaXaDw') % NI-DAQ stream: MN,MA,XA,DW info.stream_type = 'nidq'; @@ -89,14 +101,34 @@ info.n_mn_chans = counts(1); % multiplexed neural info.n_ma_chans = counts(2); % multiplexed analog info.n_xa_chans = counts(3); % non-multiplexed analog - info.n_dw_chans = counts(4); % digital words + info.n_dw_chans = counts(4); % digital word int16 columns info.n_neural_chans = counts(1) + counts(2) + counts(3); info.n_sync_chans = counts(4); + info.n_digital_word_cols = counts(4); + % SpikeGLX packs niXDBytes1 bytes from port0 plus niXDBytes2 bytes + % from port1 into the DW int16 columns. The number of meaningful + % digital lines is 8 * (niXDBytes1 + niXDBytes2). When those + % fields are absent fall back to assuming every bit of every DW + % column is in use. + n_xd_bytes = 0; + if isfield(meta, 'niXDBytes1') + n_xd_bytes = n_xd_bytes + str2double(meta.niXDBytes1); + end + if isfield(meta, 'niXDBytes2') + n_xd_bytes = n_xd_bytes + str2double(meta.niXDBytes2); + end + if n_xd_bytes > 0 + info.n_digital_lines = 8 * n_xd_bytes; + else + info.n_digital_lines = 16 * info.n_dw_chans; + end else % Fallback info.stream_type = 'unknown'; info.n_neural_chans = info.n_saved_chans - 1; info.n_sync_chans = 1; + info.n_digital_word_cols = 1; + info.n_digital_lines = 16; end % Parse saved channel subset diff --git a/+ndr/+reader/neuropixelsGLX.m b/+ndr/+reader/neuropixelsGLX.m index dfef9c5..2de3a30 100644 --- a/+ndr/+reader/neuropixelsGLX.m +++ b/+ndr/+reader/neuropixelsGLX.m @@ -12,7 +12,11 @@ % % Channel mapping: % - Neural channels are exposed as 'analog_in' (ai1..aiN) -% - The sync word is exposed as 'digital_in' (di1) +% - Digital lines are exposed as 'digital_in' (di1..diM), where each +% di channel is a single bit of the packed digital word(s). The +% number of lines is determined from metadata: for NIDQ streams +% it is 8 * (niXDBytes1 + niXDBytes2); for IMEC streams it is +% 16 * n_sync_chans (bit 6 is the SMA sync input in practice). % - A single time channel 't1' is always present % % Data is returned as int16 to preserve native precision. Use @@ -90,8 +94,9 @@ % % Returns a structure array with fields 'name', 'type', and % 'time_channel'. Neural channels are 'analog_in' (ai1..aiN), - % the sync channel is 'digital_in' (di1), and a time channel - % 't1' is always present. + % digital lines are 'digital_in' (di1..diM) with one entry + % per single-bit line in the packed digital word(s), and a + % time channel 't1' is always present. % % See also: ndr.format.neuropixelsGLX.header @@ -109,10 +114,12 @@ 'type', 'analog_in', 'time_channel', 1); %#ok end - % Sync channel (digital_in) - if info.n_sync_chans > 0 - channels(end+1) = struct('name', 'di1', ... - 'type', 'digital_in', 'time_channel', 1); + % Digital lines (digital_in) — one per bit of the packed + % digital word(s). n_digital_lines comes from metadata + % (niXDBytes1/niXDBytes2 for NIDQ, 16*n_sync_chans for IMEC). + for i = 1:info.n_digital_lines + channels(end+1) = struct('name', ['di' int2str(i)], ... + 'type', 'digital_in', 'time_channel', 1); %#ok end end @@ -167,7 +174,9 @@ % % For 'analog_in': returns int16 neural data. % For 'time': returns double time stamps in seconds. - % For 'digital_in': returns int16 sync word values. + % For 'digital_in': returns int16 single-bit values (0 or 1) + % extracted from the packed digital word(s). + % CHANNEL gives the 1-based digital line(s). % % See also: ndr.format.neuropixelsGLX.read @@ -189,9 +198,33 @@ data = read_samples(binfile, info, uint32(channel), s0, s1); case {'digital_in', 'di'} - % Sync channel is the last channel in the file - sync_chan = info.n_saved_chans; - data = read_samples(binfile, info, uint32(sync_chan), s0, s1); + % Digital words occupy the last n_digital_word_cols + % columns of the file. Each int16 column holds up to + % 16 single-bit lines. CHANNEL is a vector of 1-based + % digital line indices; map each to (column, bit) and + % extract the bit with bitget. + line_0based = double(channel(:)) - 1; + if any(line_0based < 0) || ... + any(line_0based >= info.n_digital_lines) + error('ndr:reader:neuropixelsGLX:DigitalLineOutOfRange', ... + 'Digital line out of range; valid lines are 1..%d.', ... + info.n_digital_lines); + end + first_dw_col = info.n_saved_chans - info.n_digital_word_cols + 1; + dw_col_offset = floor(line_0based / 16); % 0-based DW column offset + bit_pos = mod(line_0based, 16); % 0-based bit within column + + n_samples = double(s1) - double(s0) + 1; + data = zeros(n_samples, numel(channel), 'int16'); + unique_cols = unique(dw_col_offset); + for u = 1:numel(unique_cols) + file_col = first_dw_col + unique_cols(u); + raw = read_samples(binfile, info, uint32(file_col), s0, s1); + idx = find(dw_col_offset == unique_cols(u)); + for k = 1:numel(idx) + data(:, idx(k)) = int16(bitget(raw, bit_pos(idx(k)) + 1)); + end + end otherwise error('ndr:reader:neuropixelsGLX:UnknownChannelType', ... diff --git a/+ndr/+test/+reader/+neuropixelsGLX/test.m b/+ndr/+test/+reader/+neuropixelsGLX/test.m index 80f7bd3..eda3e26 100644 --- a/+ndr/+test/+reader/+neuropixelsGLX/test.m +++ b/+ndr/+test/+reader/+neuropixelsGLX/test.m @@ -42,11 +42,12 @@ function test(varargin) for c = 1:nNeuralChans data(:, c) = int16(round(500 * sin(2 * pi * c * t_vec))); end - % Sync channel is a 16-bit digital word (packed digital lines), not 0/1. - % Use a ramping pattern that exercises multiple bits to confirm the full - % 16-bit value is preserved end-to-end. - sync_data = int16(mod((0:nSamples-1), 2^15)); - data(:, nTotalChans) = sync_data(:); + % Sync channel is a 16-bit packed digital word. Each bit is exposed as + % a separate digital line (di1..di16). Use a counter pattern that + % exercises bits 0..14 (mod 2^15 keeps values non-negative so the int16 + % representation is unambiguous). + sync_pattern = int16(mod((0:nSamples-1), 2^15)); + data(:, nTotalChans) = sync_pattern(:); % Write binary fid = fopen(binfile, 'w', 'ieee-le'); @@ -107,23 +108,127 @@ function test(varargin) disp(['Max error vs expected: ' num2str(max_error)]); assert(max_error == 0, 'Data mismatch detected!'); - % Read the digital (sync) channel through the high-level read() API. + % Verify that the IMEC sync word is exposed as 16 single-bit digital + % lines (di1..di16), one per bit of the int16 sync column. + di_idx = find(strcmp({channels.type}, 'digital_in')); + assert(numel(di_idx) == 16, ... + sprintf('Expected 16 IMEC digital lines, got %d.', numel(di_idx))); + + % Read individual digital lines through the high-level read() API. % This is the code path that previously returned [] because % ndr.reader.read() routed 'digital_in' to readevents_epochsamples % (which is abstract for this format). - t0 = 0; - t1 = (nSamples-1) / SR; - [d_di, t_di] = r.read({metafile}, 'di1', 't0', t0, 't1', t1); - disp(['Read ' int2str(size(d_di, 1)) ' samples from channel di1 via r.read().']); - assert(~isempty(d_di), 'Digital read returned empty data.'); - assert(~isempty(t_di), 'Digital read returned empty time.'); - assert(size(d_di, 1) == nSamples, ... + t1_end = (nSamples-1) / SR; + [d_di1, t_di] = r.read({metafile}, 'di1', 't0', 0, 't1', t1_end); + disp(['Read ' int2str(size(d_di1, 1)) ' samples from channel di1 via r.read().']); + assert(~isempty(d_di1), 'Digital di1 read returned empty data.'); + assert(~isempty(t_di), 'Digital di1 read returned empty time.'); + assert(size(d_di1, 1) == nSamples, ... sprintf('Digital sample count mismatch: got %d, expected %d.', ... - size(d_di, 1), nSamples)); - % The full 16-bit word must be preserved (not collapsed to 0/1). - di_error = max(abs(double(d_di) - double(sync_data(:)))); - disp(['Max digital error vs expected: ' num2str(di_error)]); - assert(di_error == 0, 'Digital (16-bit) word mismatch detected!'); + size(d_di1, 1), nSamples)); + % di1 is bit 0 of the sync word (alternates 0,1,0,1,...). + expected_di1 = int16(bitget(sync_pattern(:), 1)); + assert(isequal(d_di1, expected_di1), 'di1 (bit 0) extraction mismatch.'); + + % di8 is bit 7 of the sync word. + [d_di8, ~] = r.read({metafile}, 'di8', 't0', 0, 't1', t1_end); + expected_di8 = int16(bitget(sync_pattern(:), 8)); + assert(isequal(d_di8, expected_di8), 'di8 (bit 7) extraction mismatch.'); + + % di12 is bit 11 of the sync word (2^11 = 2048; the pattern reaches + % 2999 so this bit is non-trivially exercised). + [d_di12, ~] = r.read({metafile}, 'di12', 't0', 0, 't1', t1_end); + expected_di12 = int16(bitget(sync_pattern(:), 12)); + assert(any(expected_di12 ~= 0), ... + 'Test pattern should exercise bit 11; check sync_pattern.'); + assert(isequal(d_di12, expected_di12), 'di12 (bit 11) extraction mismatch.'); + + % === NIDQ format test === + % Mirrors the user-reported configuration: snsMnMaXaDw=0,0,8,1 with + % niXDBytes1=1, so 8 analog inputs plus 8 digital lines packed into + % the low byte of a single int16 DW column. + disp('--- NIDQ format test ---'); + + nXA = 8; + nDW = 1; + nNidqChans = nXA + nDW; + nNidqSamples = 500; + SR_ni = 10593.220339; + + nidq_subdir = fullfile(tempdir_path, 'nidq_g0'); + mkdir(nidq_subdir); + nidq_metafile = fullfile(nidq_subdir, 'nidq_g0_t0.nidq.meta'); + nidq_binfile = fullfile(nidq_subdir, 'nidq_g0_t0.nidq.bin'); + + % Build data: XA0..XA7 sine waves + DW column with an 8-bit pattern. + ni_data = zeros(nNidqSamples, nNidqChans, 'int16'); + t_vec_ni = (0:nNidqSamples-1)' / SR_ni; + for c = 1:nXA + ni_data(:, c) = int16(round(1000 * sin(2 * pi * c * t_vec_ni))); + end + ni_sync = int16(mod((0:nNidqSamples-1), 2^8)); % low byte: 0..255 + ni_data(:, end) = ni_sync(:); + + fid = fopen(nidq_binfile, 'w', 'ieee-le'); + fwrite(fid, reshape(ni_data', 1, []), 'int16'); + fclose(fid); + + fid = fopen(nidq_metafile, 'w'); + fprintf(fid, 'niSampRate=%.6f\n', SR_ni); + fprintf(fid, 'nSavedChans=%d\n', nNidqChans); + fprintf(fid, 'snsMnMaXaDw=0,0,%d,%d\n', nXA, nDW); + fprintf(fid, 'snsSaveChanSubset=all\n'); + fprintf(fid, 'fileSizeBytes=%d\n', nNidqSamples * nNidqChans * 2); + fprintf(fid, 'fileTimeSecs=%.6f\n', nNidqSamples / SR_ni); + fprintf(fid, 'niAiRangeMax=5\n'); + fprintf(fid, 'niAiRangeMin=-5\n'); + fprintf(fid, 'niMaxInt=32768\n'); + fprintf(fid, 'niXDBytes1=1\n'); + fprintf(fid, 'niXDChans1=0:7\n'); + fprintf(fid, 'niXAChans1=0:7\n'); + fprintf(fid, 'typeThis=nidq\n'); + fclose(fid); + + r_ni = ndr.reader('neuropixelsGLX'); + ni_channels = r_ni.getchannelsepoch({nidq_metafile}, 1); + ni_di_count = sum(strcmp({ni_channels.type}, 'digital_in')); + assert(ni_di_count == 8, ... + sprintf('Expected 8 NIDQ digital lines (niXDBytes1=1), got %d.', ... + ni_di_count)); + ni_ai_count = sum(strcmp({ni_channels.type}, 'analog_in')); + assert(ni_ai_count == nXA, ... + sprintf('Expected %d NIDQ analog lines, got %d.', nXA, ni_ai_count)); + + % Read di1 via r.read() — the exact call the user reported failing. + t1_ni = (nNidqSamples-1) / SR_ni; + [d_ni_di1, t_ni_di1] = r_ni.read({nidq_metafile}, 'di1', 't0', 0, 't1', t1_ni); + assert(~isempty(d_ni_di1), 'NIDQ di1 returned empty data (user-reported bug).'); + assert(~isempty(t_ni_di1), 'NIDQ di1 returned empty time.'); + assert(size(d_ni_di1, 1) == nNidqSamples, ... + sprintf('NIDQ di1 sample count mismatch: got %d, expected %d.', ... + size(d_ni_di1, 1), nNidqSamples)); + expected_ni_di1 = int16(bitget(ni_sync(:), 1)); + assert(isequal(d_ni_di1, expected_ni_di1), 'NIDQ di1 (bit 0) mismatch.'); + + % Read di8 (bit 7, the high bit of the byte). + [d_ni_di8, ~] = r_ni.read({nidq_metafile}, 'di8', 't0', 0, 't1', t1_ni); + expected_ni_di8 = int16(bitget(ni_sync(:), 8)); + assert(isequal(d_ni_di8, expected_ni_di8), 'NIDQ di8 (bit 7) mismatch.'); + + % di9 must be out of range for an 8-line NIDQ configuration. Call + % readchannels_epochsamples directly because the high-level read() + % would fail earlier at channel-name lookup (di9 isn't listed). + threw = false; + try + r_ni.readchannels_epochsamples('digital_in', 9, ... + {nidq_metafile}, 1, 1, 10); + catch ME + threw = strcmp(ME.identifier, ... + 'ndr:reader:neuropixelsGLX:DigitalLineOutOfRange'); + end + assert(threw, 'Expected DigitalLineOutOfRange error for line 9 in 8-line NIDQ.'); + + disp('NIDQ digital line read: OK.'); disp('All checks passed.'); diff --git a/tools/tests/+ndr/+unittest/+reader/TestNeuropixelsGLX.m b/tools/tests/+ndr/+unittest/+reader/TestNeuropixelsGLX.m index 4bef82f..dddec10 100644 --- a/tools/tests/+ndr/+unittest/+reader/TestNeuropixelsGLX.m +++ b/tools/tests/+ndr/+unittest/+reader/TestNeuropixelsGLX.m @@ -160,8 +160,10 @@ function testGetChannelsEpoch(testCase) %TESTGETCHANNELSEPOCH Verify channel listing. channels = testCase.Reader.getchannelsepoch({testCase.MetaFilename}, 1); - % Should have: 1 time + N neural + 1 sync = N+2 - expectedTotal = testCase.NumNeuralChansActual + 2; + % IMEC sync word is a single int16 column = 16 single-bit + % digital lines (di1..di16). Total: 1 time + N neural + 16 di. + nDigitalLines = 16; + expectedTotal = testCase.NumNeuralChansActual + 1 + nDigitalLines; testCase.verifyNumElements(channels, expectedTotal, 'Wrong number of channels.'); % First channel should be time @@ -176,9 +178,15 @@ function testGetChannelsEpoch(testCase) ['Wrong type for neural channel ' int2str(i)]); end - % Sync channel (last) - testCase.verifyEqual(channels(end).name, 'di1', 'Last channel should be di1.'); - testCase.verifyEqual(channels(end).type, 'digital_in', 'Last channel type should be digital_in.'); + % Digital lines: di1..di16 follow the neural channels. + di_start = 1 + testCase.NumNeuralChansActual + 1; % after time + neural + for i = 1:nDigitalLines + idx = di_start + i - 1; + testCase.verifyEqual(channels(idx).name, ['di' int2str(i)], ... + ['Wrong name for digital line ' int2str(i)]); + testCase.verifyEqual(channels(idx).type, 'digital_in', ... + ['Wrong type for digital line ' int2str(i)]); + end end function testSampleRate(testCase) From c17191a77ab206bdda77885725dbc1302c94fcde Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 10 Apr 2026 15:02:26 +0000 Subject: [PATCH 3/3] Pre-compute digital line (col, bit) map and label in header The neuropixelsGLX header now stores explicit per-line metadata for the digital channels it exposes: - n_digital_lines : count of active single-bit lines - digital_line_col : 0-based DW column offset for each line - digital_line_bit : 0-based bit position within that column - digital_line_label : underlying SpikeGLX line label, e.g. 'XD0'..'XD7' for port-0, 'XD1.' for port-1, 'SY.' for IMEC sync. NIDQ line counts come from niXDBytes1/niXDBytes2 (8 active lines per saved byte). NI hardware only enables digital input in whole-byte chunks, so every bit within a captured byte is electrically active even when the user only physically wired some of them. niXDChans1/2 remains informational and is not used to gate which lines are exposed. IMEC streams have no per-bit configuration, so all 16 bits of each sync int16 column are exposed. readchannels_epochsamples no longer recomputes (col, bit) on the fly; it just looks up info.digital_line_col / digital_line_bit by the requested 1-based di index. The bit-extraction loop, the unique column grouping, and the read_samples call sites are unchanged, so behavior for the user's file (8 di channels mapping to bits 0..7 of the only DW column) is identical. --- +ndr/+format/+neuropixelsGLX/header.m | 97 ++++++++++++++++++++++----- +ndr/+reader/neuropixelsGLX.m | 25 +++---- 2 files changed, 94 insertions(+), 28 deletions(-) diff --git a/+ndr/+format/+neuropixelsGLX/header.m b/+ndr/+format/+neuropixelsGLX/header.m index 0600fde..f2b3299 100644 --- a/+ndr/+format/+neuropixelsGLX/header.m +++ b/+ndr/+format/+neuropixelsGLX/header.m @@ -68,12 +68,26 @@ info.n_saved_chans = str2double(meta.nSavedChans); % Parse snsApLfSy or snsMnMaXaDw to determine neural vs sync channels. - % Also compute: + % Also compute, for digital lines: % n_digital_word_cols : number of int16 columns in the .bin file - % that hold digital word data (stored last). + % that hold digital word data (stored last). % n_digital_lines : number of single-bit digital lines exposed. - % SpikeGLX packs digital data into int16 columns; - % each bit is an independent digital line. + % digital_line_col : (n_digital_lines x 1) 0-based DW column + % offset (0 = first DW column). + % digital_line_bit : (n_digital_lines x 1) 0-based bit position + % within that column (0..15). + % digital_line_label : (n_digital_lines x 1) cellstr describing + % the underlying SpikeGLX line, e.g. 'XD0' + % for port-0 line 0, 'XD1.3' for port-1 line + % 3, or 'SY0.6' for sync col 0 bit 6. + % + % For NIDQ streams the count of active lines comes from niXDBytes1 + % and niXDBytes2 (bytes captured per port). NI-DAQ hardware only + % enables digital input in whole-byte chunks, so every bit within a + % captured byte is electrically active even if the user only wired + % some of them; niXDChans1/2 is just informational and is not used + % to gate which lines are exposed. For IMEC streams there is no + % per-bit configuration; all 16 bits of each sync int16 are exposed. if isfield(meta, 'snsApLfSy') % imec stream: AP,LF,SY counts counts = sscanf(meta.snsApLfSy, '%d,%d,%d'); @@ -93,7 +107,20 @@ % are exposed as independent digital lines so callers can pick % whichever they need. info.n_digital_word_cols = info.n_sync_chans; - info.n_digital_lines = 16 * info.n_sync_chans; + n_lines = 16 * info.n_sync_chans; + info.n_digital_lines = n_lines; + info.digital_line_col = zeros(n_lines, 1); + info.digital_line_bit = zeros(n_lines, 1); + info.digital_line_label = cell(n_lines, 1); + idx = 0; + for c = 0:(info.n_sync_chans - 1) + for b = 0:15 + idx = idx + 1; + info.digital_line_col(idx) = c; + info.digital_line_bit(idx) = b; + info.digital_line_label{idx} = sprintf('SY%d.%d', c, b); + end + end elseif isfield(meta, 'snsMnMaXaDw') % NI-DAQ stream: MN,MA,XA,DW info.stream_type = 'nidq'; @@ -105,22 +132,54 @@ info.n_neural_chans = counts(1) + counts(2) + counts(3); info.n_sync_chans = counts(4); info.n_digital_word_cols = counts(4); - % SpikeGLX packs niXDBytes1 bytes from port0 plus niXDBytes2 bytes - % from port1 into the DW int16 columns. The number of meaningful - % digital lines is 8 * (niXDBytes1 + niXDBytes2). When those - % fields are absent fall back to assuming every bit of every DW - % column is in use. - n_xd_bytes = 0; + + % Bytes saved per port. NI hardware only enables digital input + % in whole-byte chunks, so each saved byte = 8 active lines + % regardless of how many of them the user actually wired up. + n_bytes_p0 = 0; if isfield(meta, 'niXDBytes1') - n_xd_bytes = n_xd_bytes + str2double(meta.niXDBytes1); + n_bytes_p0 = str2double(meta.niXDBytes1); end + n_bytes_p1 = 0; if isfield(meta, 'niXDBytes2') - n_xd_bytes = n_xd_bytes + str2double(meta.niXDBytes2); + n_bytes_p1 = str2double(meta.niXDBytes2); end - if n_xd_bytes > 0 - info.n_digital_lines = 8 * n_xd_bytes; + + % If neither byte field is present, fall back to assuming every + % bit of every DW int16 column is in use (16 lines per column). + if n_bytes_p0 == 0 && n_bytes_p1 == 0 + n_lines_p0 = 16 * info.n_dw_chans; + n_lines_p1 = 0; else - info.n_digital_lines = 16 * info.n_dw_chans; + n_lines_p0 = 8 * n_bytes_p0; + n_lines_p1 = 8 * n_bytes_p1; + end + + % Compute the (col, bit) position of each active line. + % SpikeGLX storage layout: port0 lines occupy the first + % n_bytes_p0*8 bits of the concatenated digital bit stream, + % then port1 lines occupy the next n_bytes_p1*8 bits. The bit + % stream is laid out across n_dw_chans int16 columns (16 bits + % per column). + n_lines = n_lines_p0 + n_lines_p1; + info.n_digital_lines = n_lines; + info.digital_line_col = zeros(n_lines, 1); + info.digital_line_bit = zeros(n_lines, 1); + info.digital_line_label = cell(n_lines, 1); + idx = 0; + for k = 0:(n_lines_p0 - 1) + abs_bit = k; + idx = idx + 1; + info.digital_line_col(idx) = floor(abs_bit / 16); + info.digital_line_bit(idx) = mod(abs_bit, 16); + info.digital_line_label{idx} = sprintf('XD%d', k); + end + for k = 0:(n_lines_p1 - 1) + abs_bit = n_bytes_p0 * 8 + k; + idx = idx + 1; + info.digital_line_col(idx) = floor(abs_bit / 16); + info.digital_line_bit(idx) = mod(abs_bit, 16); + info.digital_line_label{idx} = sprintf('XD1.%d', k); end else % Fallback @@ -129,6 +188,12 @@ info.n_sync_chans = 1; info.n_digital_word_cols = 1; info.n_digital_lines = 16; + info.digital_line_col = zeros(16, 1); + info.digital_line_bit = (0:15)'; + info.digital_line_label = cell(16, 1); + for b = 0:15 + info.digital_line_label{b+1} = sprintf('bit%d', b); + end end % Parse saved channel subset diff --git a/+ndr/+reader/neuropixelsGLX.m b/+ndr/+reader/neuropixelsGLX.m index 2de3a30..c6a01df 100644 --- a/+ndr/+reader/neuropixelsGLX.m +++ b/+ndr/+reader/neuropixelsGLX.m @@ -198,29 +198,30 @@ data = read_samples(binfile, info, uint32(channel), s0, s1); case {'digital_in', 'di'} - % Digital words occupy the last n_digital_word_cols - % columns of the file. Each int16 column holds up to - % 16 single-bit lines. CHANNEL is a vector of 1-based - % digital line indices; map each to (column, bit) and - % extract the bit with bitget. - line_0based = double(channel(:)) - 1; - if any(line_0based < 0) || ... - any(line_0based >= info.n_digital_lines) + % CHANNEL is a vector of 1-based digital line indices + % (di1..di_n_digital_lines). The header pre-computes + % the (DW column, bit position) for each active line + % from niXDBytes1/niXDBytes2 (NIDQ) or n_sync_chans + % (IMEC), so the reader just looks up the mapping + % and extracts the requested bits with bitget. + line_idx = double(channel(:)); + if any(line_idx < 1) || ... + any(line_idx > info.n_digital_lines) error('ndr:reader:neuropixelsGLX:DigitalLineOutOfRange', ... 'Digital line out of range; valid lines are 1..%d.', ... info.n_digital_lines); end first_dw_col = info.n_saved_chans - info.n_digital_word_cols + 1; - dw_col_offset = floor(line_0based / 16); % 0-based DW column offset - bit_pos = mod(line_0based, 16); % 0-based bit within column + col_offsets = info.digital_line_col(line_idx); + bit_pos = info.digital_line_bit(line_idx); n_samples = double(s1) - double(s0) + 1; data = zeros(n_samples, numel(channel), 'int16'); - unique_cols = unique(dw_col_offset); + unique_cols = unique(col_offsets); for u = 1:numel(unique_cols) file_col = first_dw_col + unique_cols(u); raw = read_samples(binfile, info, uint32(file_col), s0, s1); - idx = find(dw_col_offset == unique_cols(u)); + idx = find(col_offsets == unique_cols(u)); for k = 1:numel(idx) data(:, idx(k)) = int16(bitget(raw, bit_pos(idx(k)) + 1)); end