EvaluationMotifRetentionRate.py

import music21
import matplotlib.pyplot as plt
import itertools
from collections import defaultdict, Counter
import numpy as np
import matplotlib.cm as cm
from matplotlib.colors import LogNorm
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import colormaps

def load_musicxml(file_path):
    return music21.converter.parse(file_path)

def extract_pitches(score, part_index):
    part = score.parts[part_index]
    return [note.pitch for note in part.recurse().notes if note.isNote]

def compute_pitch_diffs(pitches):
    return [pitches[i + 1].midi - pitches[i].midi for i in range(len(pitches) - 1)]

def expand_pitch_diffs(pitch_diffs, multipliers):
    return [[diff * m for diff in pitch_diffs] for m in multipliers]

def has_matching_signs(*seqs):
    for seq1, seq2 in itertools.combinations(seqs, 2):
        if not all((x > 0 and y > 0) or (x < 0 and y < 0) for x, y in zip(seq1, seq2)):
            return False
    return True

def within_tolerance(seq1, seq2, tolerance_levels):
    for tolerance, max_exceptions in tolerance_levels:
        exceptions = 0
        for x, y in zip(seq1, seq2):
            ratio = y / x if x != 0 else 0
            if not (1 / tolerance <= abs(ratio) <= tolerance):
                exceptions += 1
                if exceptions > max_exceptions:
                    break
        else:
            return tolerance
    return None


def compute_note_to_position():
    note_to_position = {}
    for pitch in music21.scale.ChromaticScale('C1').getPitches('C1', 'C8'):
        # Add the pitch name with octave to the dictionary
        nameWithOctave = pitch.nameWithOctave
        note_to_position[nameWithOctave] = pitch.midi

        # Find the enharmonic equivalents and add them to the dictionary
        enharmonics = pitch.getAllCommonEnharmonics()
        for enharmonic in enharmonics:
            enharmonic_nameWithOctave = enharmonic.nameWithOctave
            note_to_position[enharmonic_nameWithOctave] = enharmonic.midi

        # Add the base pitch name with the current octave
        base_nameWithOctave = pitch.name + str(pitch.octave)
        note_to_position[base_nameWithOctave] = pitch.midi

    return note_to_position



def analyze_and_merge_themes(themes):
    # Group themes by pitch differences
    grouped_themes = defaultdict(list)
    for theme, info in themes.items():
        grouped_themes[theme].extend(info['indices'])

    # Analyze the density of indices for each theme group
    merged_themes = {}
    for theme, indices in grouped_themes.items():
        count_indices = Counter(indices)
        most_common_indices = count_indices.most_common()  # List of (index, count) tuples
        density_info = {'indices': most_common_indices, 'total_count': len(indices)}
        merged_themes[theme] = density_info

    return merged_themes

# def analyze_index_density(themes):
#     index_counter = Counter()
#     for theme, info in themes.items():
#         index_counter.update(info['indices'])
#     return index_counter.most_common()

def detect_recurring_themes(tracks_with_pitches, voice_names, window_size=8,
                            tolerance_levels=[(1, 3), (2, 3), (3, 0)]):
    multipliers = [-3, -2, -1, -0.5, 0.5, 1, 2, 3]
    themes = defaultdict(
        lambda: {'count': 0, 'indices': [], 'match_levels': defaultdict(int), 'first_appearance': None})
    counts_per_time_index = defaultdict(lambda: {multiplier: 0 for multiplier in multipliers})
    counts_list = []


    for track_index, (pitches, pitch_diffs) in enumerate(tracks_with_pitches):
        for i in range(len(pitch_diffs) - window_size + 1):
            window = tuple(pitch_diffs[i:i + window_size])
            expanded_windows = expand_pitch_diffs(window, multipliers)
            for other_track_index, (_, other_pitch_diffs) in enumerate(tracks_with_pitches):
                if track_index != other_track_index:
                    other_window = tuple(other_pitch_diffs[i:i + window_size])
                    for multiplier, exp_window in zip(multipliers, expanded_windows):
                        match_level = within_tolerance(exp_window, other_window, tolerance_levels)
                        if match_level is not None:
                            themes[window]['count'] += 1
                            themes[window]['indices'].append(i)
                            themes[window]['match_levels'][multiplier] += 1
                            counts_per_time_index[i][multiplier] += 1
                            if not themes[window]['first_appearance']:
                                themes[window]['first_appearance'] = (track_index, i)

    for time_index in sorted(counts_per_time_index.keys()):
        print(f"Time Index {time_index}:")
        multiplier_counts = [counts_per_time_index[time_index][m] for m in multipliers]
        counts_list.append((time_index, multiplier_counts))
    # 打印到终端
    for time_index, multiplier_counts in counts_per_time_index.items():
        print(f"Time Index {time_index}:")
        for multiplier, count in zip(multipliers, multiplier_counts):
            print(f"  Multiplier {multiplier}: Total Count = {count}")
        print()

    return themes, counts_list


def merge_themes_by_density(themes, top_indices):
    merged_themes = defaultdict(lambda: defaultdict(lambda: {'count': 0, 'themes': []}))

    for theme, info in themes.items():
        for index in top_indices:
            if index in info['indices']:
                for i, idx in enumerate(info['indices']):
                    if idx == index:
                        merged_themes[index][theme]['count'] += 1
                        if info['first_appearance'][1] == idx:
                            merged_themes[index][theme]['themes'].append(info['first_appearance'])

    return merged_themes


def compute_absolute_pitches(themes, score):
    absolute_themes = {}
    for theme_diffs, info in themes.items():
        if 'first_appearance' in info and isinstance(info['first_appearance'], tuple):
            track_index, note_index = info['first_appearance']
            if track_index < len(score.parts):
                part = score.parts[track_index]
                notes = [n for n in part.recurse().notes if n.isNote]
                if note_index < len(notes):
                    first_pitch = notes[note_index].pitch
                    absolute_themes[theme_diffs] = {
                        'first_index': note_index,
                        'first_pitch': first_pitch,
                        'density': len(info['indices'])
                    }
    return absolute_themes


def revert_to_notes(pitches, pitch_diffs, index):
    sequence = [pitches[index]]
    for diff in pitch_diffs[:7]:
        next_pitch = music21.pitch.Pitch()
        next_pitch.midi = sequence[-1].midi + diff
        sequence.append(next_pitch)
    return ' '.join(p.nameWithOctave for p in sequence)

def prepare_data_for_visualization(themes):
    data = []
    for theme, info in themes.items():
        for index in info['indices']:
            data.append((index, theme, info['count']))
    return data

def compute_relative_contours(themes):
    relative_contours = {}
    for theme, info in themes.items():
        # Only consider the pitch differences (relative contour) for the theme
        diff_sequence = theme
        # Accumulate the counts of each occurrence
        density = sum([other_info['count'] for other_info in info])
        # Use the first occurrence index for the x-axis
        first_index = min(info, key=lambda x: x['first_appearance'][1])['first_appearance'][1]
        # Use the first pitch of the first occurrence for the y-axis
        first_pitch = info[0]['first_appearance'][0]
        relative_contours[diff_sequence] = {
            'density': density,
            'first_index': first_index,
            'first_pitch': first_pitch
        }
    return relative_contours

def revert_to_absolute_pitches(pitch_diffs, starting_pitch):
    """
    Convert a sequence of pitch differences back to absolute pitches.
    """
    absolute_pitches = [starting_pitch.midi]
    for diff in pitch_diffs:
        absolute_pitches.append(absolute_pitches[-1] + diff)
    return absolute_pitches


def plot_pitch_density(absolute_themes):
    # Calculate the ranges for time indices and pitches
    max_time_index = max(info['first_index'] + len(theme) for theme, info in absolute_themes.items())
    pitches = set()
    for theme, info in absolute_themes.items():
        start_pitch = info['first_pitch'].midi if isinstance(info['first_pitch'], music21.pitch.Pitch) else info['first_pitch']
        for diff in theme:
            start_pitch += diff
            pitches.add(start_pitch)
    min_pitch, max_pitch = min(pitches), max(pitches)

    time_indices = np.arange(0, max_time_index)
    pitch_indices = np.arange(min_pitch, max_pitch + 1)
    X, Y = np.meshgrid(time_indices, pitch_indices)
    Z = np.zeros_like(X, dtype=float)

    for theme, info in absolute_themes.items():
        start_pitch = info['first_pitch'].midi if isinstance(info['first_pitch'], music21.pitch.Pitch) else info['first_pitch']
        for i, diff in enumerate(theme):
            pitch_height = start_pitch + diff
            time_index = info['first_index'] + i
            if time_index < Z.shape[1] and min_pitch <= pitch_height <= max_pitch:
                Z[pitch_height - min_pitch, time_index] = max(Z[pitch_height - min_pitch, time_index], info['density'])

    # Normalize Z by the maximum value to keep it between 0 and 1
    Z_max = np.max(Z)
    if Z_max > 0:
        Z /= Z_max

    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    surf = ax.plot_surface(X, Y, Z, cmap=cm.viridis, edgecolor='none', alpha=0.7)
    ax.set_xlabel('Time Index')
    ax.set_ylabel('Pitch')
    ax.set_zlabel('Density')
    plt.show()


def plot_time_index_vs_counts(counts_list, multipliers):
    plt.figure(figsize=(15, 7))

    # Iterate through each multiplier
    for multiplier_index, multiplier in enumerate(multipliers):
        time_indices = [item[0] for item in counts_list]
        counts = [item[1][multiplier_index] for item in counts_list]

        # Plot counts for each multiplier
        plt.plot(time_indices, counts, marker='', linestyle='-', label=f'Multiplier {multiplier}')

    plt.xlabel("Time Index")
    plt.ylabel("Total Count")
    plt.title("Total Count for Each Multiplier Over Time")
    plt.legend()
    plt.grid(True)  # Optional: to add grid lines for better readability
    plt.show()


def plot_pitch_density_2d(absolute_themes, score):
    # Create a mapping from pitch to y-axis position
    note_to_position = compute_note_to_position()

    # Initialize a matrix for pitch occurrence density
    max_time_index = max(info['first_index'] + len(theme) for theme, info in absolute_themes.items())
    density_matrix = np.zeros((len(note_to_position), max_time_index))

    # Store data for average curve calculation
    avg_curve_data = defaultdict(lambda: {'total_weight': 0, 'weighted_pitches': np.zeros(8)})

    # Populate the density matrix and calculate average curves
    for theme, info in absolute_themes.items():
        start_pitch = info['first_pitch']
        start_position = note_to_position[start_pitch.nameWithOctave]
        density = info['density']
        time_index = info['first_index']

        # Calculate average curve data
        for i, diff in enumerate(theme[:8]):
            next_pitch = start_pitch.transpose(diff)
            next_position = note_to_position[next_pitch.nameWithOctave]
            if 0 <= time_index + i < max_time_index:
                density_matrix[next_position, time_index + i] += 1
                avg_curve_data[time_index]['total_weight'] += density
                avg_curve_data[time_index]['weighted_pitches'][i] += density * next_position

    # Normalize average curves
    for time_index, data in avg_curve_data.items():
        if data['total_weight'] > 0:
            data['weighted_pitches'] /= data['total_weight']

    # Plotting
    plt.figure(figsize=(12, 6))
    ax = plt.gca()
    ax.set_facecolor('black')

    # Plot average curves
    for time_index, data in avg_curve_data.items():
        x_vals = np.arange(time_index, time_index + 8)
        y_vals = data['weighted_pitches']
        ax.plot(x_vals, y_vals, linewidth=2, color='cyan')

    # Plot the density matrix
    density_log_scaled = np.log(density_matrix + 1)
    plt.imshow(density_log_scaled, aspect='auto', origin='lower', cmap='Greens',
               extent=[0, max_time_index, 35, 85],
               norm=LogNorm(vmin=np.min(density_log_scaled[np.nonzero(density_log_scaled)]),
                            vmax=np.max(density_log_scaled)))
    plt.colorbar(label='Log-scaled Density')
    plt.xlabel('Time Index')
    plt.ylabel('Pitch')
    plt.title('Average Melody Patterns by Density (2D Scroll-Like Plot)')
    plt.show()


def convert_midi_to_score(pitch_midi):
    """
    Convert MIDI pitch numbers to score notation.
    """
    return music21.pitch.Pitch(midi=pitch_midi).nameWithOctave
def weighted_average_pitch(pitches, weights):
    return np.average(pitches, weights=weights)

def plot_grouped_motifs(absolute_themes, score):
    # Compute the mapping from pitch to y-axis position
    note_to_position = compute_note_to_position()

    # Determine the range for time indices
    max_time_index = max(info['first_index'] + len(theme) for theme, info in absolute_themes.items())
    time_indices = np.arange(0, max_time_index)

    # Initialize data structure for average curve calculation
    avg_curve_data = defaultdict(lambda: {'weights': [], 'pitches': []})

    # Collect data for average curve calculation
    for theme, info in absolute_themes.items():
        start_pitch = info['first_pitch']
        density = info['density']
        time_index = info['first_index']

        # Only consider the first 8 pitch differences
        for i, diff in enumerate(theme[:8]):
            next_pitch = start_pitch.transpose(diff).midi
            position = note_to_position[start_pitch.nameWithOctave] + i
            avg_curve_data[time_index]['weights'].append(density)
            avg_curve_data[time_index]['pitches'].append(next_pitch)

    # Plot the grouped lines
    plt.figure(figsize=(12, 6))
    for time_index in range(0, max_time_index, 4):  # Plot one line per 4 time ticks
        weights = []
        pitches = []
        for offset in range(-4, 4):  # Consider an 8-time tick range
            if time_index + offset in avg_curve_data:
                weights.extend(avg_curve_data[time_index + offset]['weights'])
                pitches.extend(avg_curve_data[time_index + offset]['pitches'])

        if weights:  # If there are motifs to plot
            avg_pitch = weighted_average_pitch(pitches, weights)
            plt.plot(range(time_index, time_index + 8), [avg_pitch] * 8, label=f'Time index: {time_index}')

    plt.xlabel('Time Index')
    plt.ylabel('Pitch')
    plt.title('Grouped Motifs Over Time')
    plt.legend()
    plt.show()

def plot_weighted_motif_curves(terminal_outputs):
    # Initialize structures to hold data
    time_index_to_pitches = defaultdict(lambda: defaultdict(int))
    pitch_counts = defaultdict(Counter)

    # Parse the terminal outputs to populate time_index_to_pitches
    for output in terminal_outputs:
        initial_index = output['Initial Index']
        count = output['Count']
        pitch_diffs = output['Pitch Difference']
        theme = output['Theme']

        # Assume theme is a list of pitches; convert to MIDI if necessary
        pitches = [music21.pitch.Pitch(p) for p in theme.split()]
        for i, pitch in enumerate(pitches):
            time_index_to_pitches[initial_index + i][pitch.midi] += count

    # Calculate the weighted average pitch for each time index
    for time_index, pitches in time_index_to_pitches.items():
        for pitch, count in pitches.items():
            pitch_counts[time_index].update({pitch: count})

    # Determine the most common time index (t_highest) to start plotting from
    t_highest_counts = Counter({idx: sum(pitch_counts[idx].values()) for idx in pitch_counts})
    t_highest = t_highest_counts.most_common(1)[0][0]

    # Prepare the plot
    plt.figure(figsize=(12, 6))

    # Plot up to three lines starting from t_highest to t_highest + 8
    for start_idx in range(t_highest, t_highest + 8, 4):  # Adjust step if necessary
        for offset in range(3):  # Maximum of 3 lines
            if start_idx + offset in pitch_counts:
                time_range = range(start_idx + offset, start_idx + offset + 8)
                weighted_pitches = [
                    np.average(
                        list(pitch_counts[start_idx + offset].elements()),
                        weights=list(pitch_counts[start_idx + offset].values())
                    ) for _ in time_range
                ]
                plt.plot(time_range, weighted_pitches, label=f'Time index: {start_idx + offset}')

    # Finalize the plot
    plt.xlabel('Time Index')
    plt.ylabel('Weighted Average Pitch')
    plt.title('Weighted Average Pitch Curves for Motifs Over Time')
    plt.legend()
    plt.show()

# Example usage:
# Assuming `absolute_themes` and `score` are already defined
#plot_melodies_on_scroll(absolute_themes, score, top_n=10)
def main(file_path):
    score = load_musicxml(file_path)
    voice_names = ['Soprano', 'Alto', 'Tenor', 'Bass']
    cross_voice_ranges = [('Soprano-Alto Cross', (65, 70.5)), ('Alto-Tenor Cross', (58.5, 65)), ('Tenor-Bass Cross', (50.5, 58.5))]
    multipliers = [-3, -2, -1, -0.5, 0.5, 1, 2, 3]
    all_tracks_with_pitches = []
    for index, voice_name in enumerate(voice_names):
        pitches = extract_pitches(score, index)
        all_tracks_with_pitches.append((pitches, compute_pitch_diffs(pitches)))

    for cross_voice_name, (lower, upper) in cross_voice_ranges:
        cross_voice_pitches = []
        for part in score.parts:
            cross_voice_pitches.extend([note.pitch for note in part.recurse().notes if lower <= note.pitch.midi <= upper])
        all_tracks_with_pitches.append((cross_voice_pitches, compute_pitch_diffs(cross_voice_pitches)))

        # 正确的调用和赋值
    themes, counts_list = detect_recurring_themes(all_tracks_with_pitches, voice_names)

    # 打印每个时间索引的计数信息
    print("计数每个时间索引:")
    for time_index, multiplier_counts in counts_list:
        print(f"时间索引 {time_index}:")
        for multiplier_index, count in enumerate(multiplier_counts):
            print(f"  Multiplier {multipliers[multiplier_index]}: Total Count = {count}")
        print()

    # After detecting recurring themes:
    theme_data = detect_recurring_themes(all_tracks_with_pitches, voice_names)

    themes = theme_data[0]  # Get the themes dictionary
    counts_per_time_index = theme_data[1]  # Get the counts_per_time_index dictionary

    # 绘制图表
    absolute_themes = compute_absolute_pitches(themes, score)
    plot_pitch_density(absolute_themes)
    plot_pitch_density_2d(absolute_themes, score)
    plot_time_index_vs_counts(counts_list, multipliers)

if __name__ == "__main__":
    file_path = r"C:\Users\Tyan\Desktop\Experiment Data\14.mxl"
    # file_path = r"C:\TrainedDataForResearch\16barGenerationGood\6.mxl"modified
    #comparison pair 1: original2 : modifed 6
    #                   origina11  :  modifed 2
    #                   original4 : modified 11
    #                   original 22: modified 14 used as demo in thesis defense




    main(file_path)