Lecture 4: Attention Detection and Music Emotion Recognition

def receive_data_in_chunks():
    print(f"Looking for an EEG stream...")
    streams = resolve_byprop("type", "EEG")

    # create a new inlet to read from stream
    inlet = StreamInlet(streams[0])
    """
    TODO:
    1. Use inlet.pull_chunk() function to receive the EEG data. 
        You can return two lists (samples and timestamps) or use parameter `dest_obj` to put EEG data directly into the buffer.
    2. Resample the data so that the data is 250 points per second using scipy.signal.resample()
    3. Print sentences to demonstrate you eeg data shape if you get the eeg data successfully.
    """
    # Get stream info to determine sampling rate and number of channels
    info = inlet.info()
    original_fs = info.nominal_srate()  # Original sampling rate
    n_channels = info.channel_count()  # Number of EEG channels
    target_fs = 250  # Target sampling rate (Hz)

    # Initialize lists to store samples and timestamps
    samples = []
    timestamps = []

    # Pull chunk of EEG data
    chunk, ts = inlet.pull_chunk(timeout=1.0, max_samples=1024)

    if chunk:
        # Convert chunk to numpy array
        samples = np.array(chunk)
        timestamps = np.array(ts)

        # Check if data was received
        if samples.size > 0:
            print(f"Received EEG data with shape: {samples.shape}")

            # Calculate number of samples for resampling
            duration = len(timestamps) / original_fs
            n_samples_new = int(duration * target_fs)

            # Resample the data to 250 Hz
            resampled_data = resample(samples, n_samples_new, axis=0)

            print(f"Resampled EEG data to 250 Hz, new shape: {resampled_data.shape}")

            return resampled_data, timestamps
        else:
            print("No EEG data received in the chunk.")
            return None, None
    else:
        print("Failed to pull EEG data chunk.")
        return None, None

    pass

def pull_and_plot(self, plot_time, plt):
    # pull the data
    _, ts = self.inlet.pull_chunk(
        timeout=0.0, max_samples=self.buffer.shape[0], dest_obj=self.buffer
    )
    # ts will be empty if no samples were pulled, a list of timestamps otherwise
    if ts:
        ts = np.asarray(ts)
        # y的形状是(样本数, 通道数)，表示多通道信号数据
        y = self.buffer[0 : ts.size, :]

        y = self.bandpass_filter(y, sfreq=250, l_freq=1.0, h_freq=45.0)
        y = self.sliding_demean_fast(y)

        # this_x：用于存储合并后的时间戳（稍后计算）。
        # old_offset和new_offset：用于确定旧数据和新数据的裁剪位置。
        this_x = None
        old_offset = 0
        new_offset = 0
        # 遍历所有通道（self.channel_count表示通道数量）
        for ch_ix in range(self.channel_count):
            # we don't pull an entire screen's worth of data, so we have to
            # trim the old data and append the new data to it
            old_x, old_y = self.curves[ch_ix].getData()
            # the timestamps are identical for all channels, so we need to do
            # this calculation only once
        if ch_ix == 0:
                # find the index of the first sample that's still visible,
                # i.e. newer than the left border of the plot
                old_offset = old_x.searchsorted(plot_time)
                # same for the new data, in case we pulled more data than
                # can be shown at once
                new_offset = ts.searchsorted(plot_time)
                # append new timestamps to the trimmed old timestamps
                this_x = np.hstack((old_x[old_offset:], ts[new_offset:]))
            # append new data to the trimmed old data
            this_y = np.hstack((old_y[old_offset:], y[new_offset:, ch_ix] - ch_ix))
            # replace the old data
            self.curves[ch_ix].setData(this_x, this_y)

def calc_bandpower(raw: mne.io.RawArray, bands=None):
    """
    Calculate the power of each channel in each frequency band
    :param raw: mne.io.RawArrray
    :param bands: Band definition dictionary, e.g. {'alpha': (8, 13), ...}
    :return: pandas.DataFrame, each row is a channel and each column is a frequency band
    """
    if bands is None:
        bands = {
            'delta': (1, 4),
            'theta': (4, 8),
            'alpha': (8, 13),
            'beta':  (13, 30),
            'gamma': (30, 45)
        }

    # compute psd
    psds_obj = raw.compute_psd(fmin=1., fmax=45., n_fft=512, verbose=False)
    psds_obj.plot(show=True)
    plt.show()
    psds = psds_obj.get_data()  # shape: (n_channels, n_freqs)
    freqs = psds_obj.freqs

    # psds.shape = (n_channels, n_freqs)
    psds = psds * 1e12  # V^2/Hz → uV^2/Hz (optional)

    """
    TODO
    """
    # DataFrame
    bandpower = {}
    for band, (fmin, fmax) in bands.items():
        freq_mask = (freqs >= fmin) & (freqs <= fmax)
        bp = psds[:, freq_mask].sum(axis=1)
        bandpower[band] = bp

    df = pd.DataFrame(bandpower, index=raw.ch_names)

    return df.round(2)