assignment4.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import scipy\n",
    "import pretty_midi as pm\n",
    "import matplotlib.pyplot as plt\n",
    "import librosa\n",
    "import librosa.display\n",
    "import math\n",
    "import itertools\n",
    "from IPython.display import Audio\n",
    "from IPython.display import Video\n",
    "import time\n",
    "import soundfile as sf\n",
    "import subprocess\n",
    "import matplotlib.animation as anim\n",
    "from moviepy.editor import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!python data_cleaning.py -i \"./data/solar_wind_data_2003-10-27 - 2003-11-02 1min.xlsx\" -o \"./data/cleaned_data.csv\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class Sonify:\n",
    "  \n",
    "  def read_data(self, input_path, output_path): \n",
    "    \"\"\"Read the data from the given csv file path.\"\"\"\n",
    "\n",
    "    try:\n",
    "      if not output_path.endswith(\".csv\"): # If provided path is not a valid csv file\n",
    "        print(\"Invalid file path. Must be .csv file.\")\n",
    "        raise FileNotFoundError()\n",
    "      else:\n",
    "        self.df = pd.read_csv(output_path) # Loads csv file\n",
    "        print(\"Dataset loaded correctly:\", output_path)\n",
    "\n",
    "      # Constants used to create Kp comparison value\n",
    "      self.density_max = np.max(self.df[\"proton_density\"])\n",
    "      self.speed_max = np.max(self.df[\"speed\"])\n",
    "      self.dataset_len = len(self.df.index)\n",
    "\n",
    "    except Exception as e:\n",
    "      if e.__class__ == FileNotFoundError:\n",
    "        print(\"No valid .csv file found at specified path. Please verify filepath.\")\n",
    "      else:\n",
    "        print(\"Error:\", e.__class__)\n",
    "\n",
    "\n",
    "  def read_midi(self, path):\n",
    "    \"\"\"Loads the MIDI file, computes grain start points as 1/32 notes.\"\"\"\n",
    "    try:\n",
    "      if not path.endswith(\".mid\"): # Checks for if the path provided is a MIDI file\n",
    "        print(\"Invalid file path. Must be .mid file.\")\n",
    "        raise FileNotFoundError()\n",
    "      else:\n",
    "        midi_data = pm.PrettyMIDI(path) # Loads the midi file\n",
    "        print('MIDI file loaded successfully:', path)\n",
    "\n",
    "      segments = (np.array(pm.PrettyMIDI.get_beats(midi_data,start_time=0.0))) # Gets the beat locations in MIDI file in seconds\n",
    "      self.midi_grain_starts = np.linspace(0, segments[-1], len(segments)*8) # Interpolates down to 1/32 note locations\n",
    "\n",
    "      # Creates attribute for MIDI grain length as the distance between steps in self.midi_grain_starts, converted to samples\n",
    "      self.midi_grain_len = int(self.midi_grain_starts[1] * self.sr)\n",
    "\n",
    "    except Exception as e:\n",
    "      if e.__class__ == FileNotFoundError:\n",
    "        print(\"MIDI file not found. Please verify filepath.\")\n",
    "      else:\n",
    "        print(\"Error:\", e.__class__)\n",
    "\n",
    "  def read_audio(self, path):\n",
    "    \"\"\"Read the audio from given filepath.\"\"\"\n",
    "    self.sr = 48000 # Initiates the sample rate for the sonification process\n",
    "\n",
    "    # Tries to load the audio file, raises FileNotFoundError if not possible.\n",
    "    try:\n",
    "      if not path.endswith(\".wav\"):\n",
    "        raise FileNotFoundError()\n",
    "      else:\n",
    "        self.song = librosa.load(path, sr=self.sr)[0]\n",
    "        print(\"Audio loaded correctly:\", path)\n",
    "\n",
    "    except Exception as e:\n",
    "      if e.__class__ == FileNotFoundError:\n",
    "        print(\"Not a valid .wav file. Please verify filepath.\")\n",
    "      else:\n",
    "        print(\"Error:\", e.__class__)\n",
    "        \n",
    "        \n",
    "\n",
    "  def __dp_grains_dataframe(self):\n",
    "    \"\"\"Creates the dataframe containing all datapoint grain IDs, start points and feature extraction values.\"\"\"\n",
    "    self.dp_grains_data = pd.DataFrame() # Data frame for grains segmented based on number of data points\n",
    "\n",
    "    song_len_samp = self.song.size # Length of current song in samples\n",
    "    self.total_grains_in_song = self.dataset_len # Number of grains in song (just the length of the dataset)\n",
    "    self.dp_grain_len = math.floor(song_len_samp / self.total_grains_in_song) # The length of each grain in samples\n",
    "\n",
    "    for i in range(self.dataset_len): # Creates one grain for each row of dataframe\n",
    "\n",
    "      self.dp_grains_data.at[i, \"grain_id\"] = i # Adds grain ID to dataframe\n",
    "\n",
    "      grain_start = i*self.dp_grain_len # Find starts sample for datapoint segmented grain\n",
    "      self.dp_grains_data.at[i, \"grain_start\"] = grain_start # Adding grain start sample index to dataframe\n",
    "      \n",
    "      # Feature extraction, returns list of feature values\n",
    "      feature_values = self.__extract_features_from_grain(grain_start, self.dp_grain_len)\n",
    "\n",
    "      # Add features to dataframe\n",
    "      self.dp_grains_data.at[i, \"spec_band\"] = feature_values[0]\n",
    "      self.dp_grains_data.at[i, \"rms\"] = feature_values[1]\n",
    "      self.dp_grains_data.at[i, \"mfcc\"] = feature_values[2]\n",
    "        \n",
    "        \n",
    "\n",
    "  def __midi_grains_dataframe(self):\n",
    "    \"\"\"Creates the dataframe containing all MIDI grain IDs, start points and feature extraction values.\"\"\"\n",
    "    self.midi_grains_data = pd.DataFrame() # Data frame for grains segmented based on previously segmented MIDI file 1/32 notes\n",
    "\n",
    "    for i in range(self.dataset_len):\n",
    "      self.midi_grains_data.at[i, \"grain_id\"] = i # Adds grain ID to dataframe\n",
    "      \n",
    "      grain_start = int(self.midi_grain_starts[i] * self.sr) # Find starts sample for MIDI segmented grain\n",
    "      self.midi_grains_data.at[i, \"grain_start\"] = grain_start # Adds grain start sample index to dataframe\n",
    "\n",
    "      # Feature extraction, returns list of feature values\n",
    "      feature_values = self.__extract_features_from_grain(grain_start, self.midi_grain_len)\n",
    "\n",
    "      # Add features to dataframe\n",
    "      self.midi_grains_data.at[i, \"spec_band\"] = feature_values[0]\n",
    "      self.midi_grains_data.at[i, \"rms\"] = feature_values[1]\n",
    "      self.midi_grains_data.at[i, \"mfcc\"] = feature_values[2]\n",
    "\n",
    "  def build_grains_dataframes(self):\n",
    "    \"\"\"Utility method to create dataframes for both types of segmentation, by datapoints and by MIDI\"\"\"\n",
    "    self.__dp_grains_dataframe()\n",
    "    self.__midi_grains_dataframe()\n",
    "\n",
    "  def __extract_features_from_grain(self, grain_start, grain_len):\n",
    "    \"\"\"Receives grain position info and computes spectral bandwidth, RMS and MFCC Librosa feature extraction. Returns results.\"\"\"\n",
    "    grain_start = int(grain_start)\n",
    "    grain_len = int(grain_len)\n",
    "    \n",
    "    grain = self.song[grain_start:grain_start+grain_len] # Selecting grain data based on start and end points\n",
    "\n",
    "    # Spectral Bandwidth\n",
    "    spec_band = librosa.feature.spectral_bandwidth(y=grain, sr=self.sr, n_fft=grain_len, hop_length=grain_len, win_length=grain_len)\n",
    "    \n",
    "    # Root Mean Square\n",
    "    rms = librosa.feature.rms(y=grain, frame_length=grain_len, hop_length=grain_len)\n",
    "    \n",
    "    # Mel Frequency Cepstrum Coefficients\n",
    "    mfcc_coeffs = librosa.feature.mfcc(y=grain, sr=self.sr, n_mels=13, n_fft=grain_len, hop_length=grain_len, win_length=grain_len)\n",
    "    mfcc = math.sqrt(np.max(mfcc_coeffs)) # Gets sqrt of largest coefficient, to make range more usable.\n",
    "\n",
    "    return [spec_band, rms, mfcc]\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def __map_feature_iterator(self, feature_list, datapoints_list, grains_dataframe, grains_type=\"dp\"):    \n",
    "    \"\"\"\"Maps the provided list of features to the provided list of datapoints in the dataset via assigning grain IDs.\n",
    "        Called inside the map_features() method.\"\"\"\n",
    "    \n",
    "    try:\n",
    "      if grains_type not in [\"dp\", \"midi\"]:\n",
    "        raise ValueError()\n",
    "\n",
    "      mapping_start_time = time.time()\n",
    "\n",
    "      for i, feature in enumerate(feature_list):\n",
    "\n",
    "        # MIDI grains data sorted by current feature\n",
    "        feature_sorted = grains_dataframe[[\"grain_id\",feature]].sort_values(by=feature, axis='index', kind=\"mergesort\", ignore_index=True)\n",
    "        # Dataset dataset by current mapped value in datapoints_list\n",
    "        self.df.sort_values(by=datapoints_list[i], axis=\"index\", kind=\"mergesort\", ignore_index=True, inplace=True)\n",
    "\n",
    "        grain_id_column_title = grains_type + \"_\" + feature + \"_grain_id\" # Assembling new grain id column title for self.df\n",
    "        feature_sorted.columns = [grain_id_column_title, feature] # Renaming columns in feature_sorted based on new titles\n",
    "\n",
    "        # Adds the new grain ID column to the dataset self.df\n",
    "        self.df = pd.concat([self.df, feature_sorted[grain_id_column_title]], axis=\"columns\")\n",
    "\n",
    "        # Re-sorts sorted self.df dataset\n",
    "        self.df.sort_values(by=\"cumulative_secs\", axis=\"index\", kind=\"mergesort\", ignore_index=True, inplace=True)\n",
    "\n",
    "      # Determines time to map features for the current grain type and prints the result\n",
    "      mapping_exec_time = time.time() - mapping_start_time\n",
    "      print(\"Mapping \", grains_type, \"grains completed in\", round(mapping_exec_time, 3), \"secs\")\n",
    "\n",
    "    except Exception as e:\n",
    "      if e.__class__ == ValueError:\n",
    "        print(\"Ensure grains_type is either 'dp' or 'midi'.\")\n",
    "      else:\n",
    "        print(\"Error:\", e.__class__)\n",
    "        \n",
    "        \n",
    "        \n",
    "  \n",
    "  def map_features(self):\n",
    "    \"\"\"Utility method which adds grain ID columns for both grain types to dataframe by calling __map_feature_iterator method\"\"\"\n",
    "    feature_list = [\"spec_band\", \"rms\", \"mfcc\"] # Features to be mapped\n",
    "    datapoints_list = [\"proton_density\", \"bz\", \"ion_temp\"] # Datapoints to be mapped\n",
    "    \n",
    "    # For datapoint grains\n",
    "    grains_type = \"dp\"\n",
    "    self.__map_feature_iterator(feature_list, datapoints_list, self.dp_grains_data, grains_type)\n",
    "\n",
    "    # For MIDI grains\n",
    "    grains_type = \"midi\"\n",
    "    self.__map_feature_iterator(feature_list, datapoints_list, self.midi_grains_data, grains_type)\n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "  def __sum_grains(self, list_of_grains):\n",
    "    \"\"\"Adds together the individual grains for each datapoint into the full grain.\"\"\"\n",
    "    # how do we want to compound the grains?\n",
    "    # just sums them for now\n",
    "    return np.sum(list_of_grains, axis=0)\n",
    "\n",
    "\n",
    "\n",
    "    \n",
    "  def scale_features(self):\n",
    "    \"\"\"Scale the bz and speed data fields from self.df dataset to required ranges as audio processing parameters.\"\"\"\n",
    "    \n",
    "    # Scales Bz to between 1 and -1\n",
    "    bz = self.df['bz'].copy()\n",
    "    bz = bz/np.max(bz)\n",
    "    self.df['bz_scaled'] = bz\n",
    "    \n",
    "    # Scales speed to between 0.5 and 2\n",
    "    speed = self.df['speed'].copy()\n",
    "    speed -= speed.min() # Scales speed to between 0 and max-min\n",
    "    speed /= speed.max() # Scales speed to between 0 and 1\n",
    "    speed = (speed * 1.5) + 0.5 # Scales speed to between 0.5 and 2\n",
    "    self.df['speed_scaled'] = speed\n",
    "    \n",
    "    print('Scaling of bz and speed values complete and written to dataframe.')\n",
    "    \n",
    "    \n",
    "  def __apply_bz(self, grain, i):\n",
    "    \"\"\"If bz is positive, apply low-pass filter. If bz is negative, apply high-pass filter.\n",
    "        Uses abs(bz) to create scaled cutoff value for use in filter.\"\"\"\n",
    "    bz_scaled = self.df.iloc[i]['bz_scaled']\n",
    "    \n",
    "    normal_cutoff = (abs(bz_scaled)*0.5) + 0.01 # Generates the cutoff value for signal.butter, in range 0.01-0.5\n",
    "    order = 5\n",
    "    b, a = scipy.signal.butter(order, normal_cutoff, btype=('low' if bz_scaled > 0 else 'high'),\n",
    "                               analog=False)\n",
    "    return scipy.signal.filtfilt(b, a, grain)\n",
    "    \n",
    "  def __apply_speed(self, grain, i, grain_len):\n",
    "    \"\"\"Change playback speed based on scaled speed value in dataset.\"\"\"\n",
    "    speed_scaled = self.df.iloc[i]['speed_scaled'] # Gets the scaled speed value.\n",
    "    grain = librosa.effects.time_stretch(grain, rate=speed_scaled) # Applies the time-stretching to the grain.\n",
    "    grain = np.append(grain, grain) # Duplicates the time-stretched grain, so it still fills the grain length if it's sped up.\n",
    "    grain = grain[:grain_len] # Trimming the duplicated grain to the grain length.\n",
    "    \n",
    "    return grain\n",
    "\n",
    "  def __apply_phi(self, grain, i):\n",
    "    \"\"\"Change playback direction based on sign of phi angle value in dataset.\"\"\"\n",
    "    phi = self.df.iloc[i]['phi_angle']\n",
    "    \n",
    "    return np.flip(grain) if phi < 0 else grain\n",
    "\n",
    "  def __simple_envelope(self, grain, grain_len):\n",
    "    \"\"\"Apply a simple envelope (100-sample-long fade on each end) to the grain.\"\"\"\n",
    "    a = np.arange(0,1,0.01) # Ramp up of the envelope\n",
    "    len_a = len(a) # The length of the ramp\n",
    "    flat = np.ones((grain_len-2 * len_a)) # Flat middle section of envelope (all 1 values)\n",
    "    rev_a = np.flip(a) # Creating ramp down of envelope by flipping the ramp up\n",
    "    env = np.append(a,flat) # Adding ramp up and flat to envelope\n",
    "    env = np.append(env,rev_a) # Adding together all the pieces of the envelope\n",
    "\n",
    "    return grain * env\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def __compound_grains(self, grains_type='dp'):\n",
    "    \"\"\"Finds, sums, processes and appends all relevant grains from the dataset.\n",
    "      Uses grains_type string to determine whether to use datapoint or MIDI segmentation.\"\"\"\n",
    "    \n",
    "    synth_start_time = time.time() # Used for measuring execution time of the loop.\n",
    "    \n",
    "    # Tries to assign the dataframe, columns and grain length based on grains_type. Generates error if invalid input.\n",
    "    try:\n",
    "      if grains_type == 'dp':\n",
    "        cols = ['dp_spec_band_grain_id', 'dp_rms_grain_id', 'dp_mfcc_grain_id']\n",
    "        grains_dataframe = self.dp_grains_data\n",
    "        grain_len = self.dp_grain_len\n",
    "      elif grains_type == 'midi':\n",
    "        cols = ['midi_spec_band_grain_id', 'midi_rms_grain_id', 'midi_mfcc_grain_id']\n",
    "        grains_dataframe = self.midi_grains_data\n",
    "        grain_len = self.midi_grain_len\n",
    "      else:\n",
    "        raise ValueError()\n",
    "\n",
    "    except ValueError:\n",
    "      print(\"grains_type must be either 'dp' or 'midi'.\")\n",
    "\n",
    "    # Main loop build the output from individual summed then processed grains\n",
    "    grains_attached = [] # Initialises main return list containing all audio\n",
    "    window = scipy.signal.windows.hamming(grain_len) # Generates the window function to be used for windowing the grain\n",
    "\n",
    "    for i, row in self.df.iterrows(): # Iterating through rows of the dataset\n",
    "      list_grains = [] # Resets list containing the data from the song to be summed to make compound grain\n",
    "      for c in cols: # Loops through columns to pull grain audio samples\n",
    "        grain_start = grains_dataframe[grains_dataframe['grain_id'] == row[c]]['grain_start'].values[0] # Gets grain based on grain ID in current column\n",
    "        grain_end = grain_start + grain_len # Calculates grain end point\n",
    "        grain = self.song[int(grain_start):int(grain_end)] # Pulls the grain audio\n",
    "        list_grains.append(grain) # Appends to list of grains for current row in data\n",
    "      compound_grain = self.__sum_grains(list_grains) # Sums the selected grains using \n",
    "        \n",
    "      # print(i, '/', self.total_grains_in_song, ' merging grains...', end='\\r')\n",
    "      \n",
    "      # On each compound grain, apply processing from the mappings\n",
    "      compound_grain = self.__apply_bz(compound_grain, i) # Applies lp/hp filter\n",
    "      compound_grain = self.__apply_speed(compound_grain, i, grain_len) # Changes playback speed\n",
    "      compound_grain = self.__apply_phi(compound_grain, i) # Changes playback direction\n",
    "      \n",
    "      # Simple envelope to ramp start and end\n",
    "      compound_grain = compound_grain * window\n",
    "      \n",
    "      # Append to main output audio array\n",
    "      grains_attached = np.append(grains_attached, compound_grain)\n",
    "    \n",
    "    synth_exec_time = time.time() - synth_start_time\n",
    "    print(\"Completed assembly of\", grains_type, \"grains in\", round(synth_exec_time, 3), \"secs\")\n",
    "\n",
    "    return grains_attached\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def synthesise_both(self):\n",
    "    \"\"\"Utility method to assemble both types of segmentation\"\"\"\n",
    "    print('--------Coffee break advised.--------')\n",
    "    \n",
    "    self.dp_synth = self.__normalise(self.__compound_grains(\"dp\"))\n",
    "    self.midi_synth = self.__normalise(self.__compound_grains(\"midi\"))\n",
    "\n",
    "    return (self.dp_synth, self.midi_synth)\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def __normalise(self, audio):\n",
    "    \"\"\"Utility method to normalise the given audio.\"\"\"\n",
    "    audio_max = np.max(np.abs(audio))\n",
    "    norm_audio = audio*1/audio_max\n",
    "    return norm_audio\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def __apply_reverb(self, path, grains_type=\"dp\"):\n",
    "    \"\"\"Applies convolution reverb audios in sections of 36 grains, as Kp-index only updates every 3 hours/36 rows of dataset.\"\"\"\n",
    "    rev_start_time = time.time()\n",
    "    \n",
    "    # Allows method to determine which audio to use based on grains_type string, raises error if not valid option\n",
    "    try:\n",
    "      if grains_type == \"dp\":\n",
    "        synth = self.dp_synth\n",
    "        grain_len = self.dp_grain_len\n",
    "      elif grains_type == \"midi\":\n",
    "        synth = self.midi_synth\n",
    "        grain_len = self.midi_grain_len\n",
    "      else:\n",
    "        raise ValueError()\n",
    "    except ValueError as e:\n",
    "      print(\"grains_type must be either 'dp' or 'midi'.\")\n",
    "      \n",
    "    output_reverb = np.zeros(synth.size)\n",
    "    kp_section_no_of_datapoints = 36 # Number of datapoints in a Kp section\n",
    "    kp_section_len = kp_section_no_of_datapoints*grain_len # Length of a Kp section in samples\n",
    "    section_window = scipy.signal.windows.tukey(kp_section_len, alpha=0.01) # Nearly rectangular Tukey window\n",
    "\n",
    "    # Tries to load the impulse response and makes sure it is shorter than one section.\n",
    "    try:\n",
    "      ir = librosa.load(path, sr=self.sr, mono=True)[0] # Loads the IR\n",
    "      ir = self.__normalise(ir) # Normalises the IR\n",
    "      if ir.size > kp_section_len: # The IR can't be longer than the kp_section, so this checks.\n",
    "        raise ValueError()\n",
    "\n",
    "    except Exception as e:\n",
    "      if e.__class__ == ValueError: # If the IR is too long, just return the un-reverbed audio\n",
    "        print(\"You must use an impulse response that is shorter than\", kp_section_len/self.sr, \".\")\n",
    "        print(\"Returning non-reverbed audio...\")\n",
    "        return synth\n",
    "      elif e.__class__ == FileNotFoundError: # If the IR isn't a valid .wav, just return the un-reverbed audio.\n",
    "        print(\"Your impulse response must be a valid .wav file.\")\n",
    "        print(\"Returning non-reverbed audio...\")\n",
    "        return synth\n",
    "      else:\n",
    "        print(\"Error:\", e.__class__)\n",
    "        print(\"Returning non-reverbed audio...\")\n",
    "        return synth\n",
    "\n",
    "    ir_padded = np.pad(ir, (0, kp_section_len - ir.size)) # Padding the IR to the same size as the section\n",
    "\n",
    "    for i in range(0, self.dataset_len, kp_section_no_of_datapoints): # Loop for number of sections\n",
    "      kp_section_start = i*grain_len # Start point in samples of current section\n",
    "      kp_section = synth[kp_section_start:kp_section_start + kp_section_len] # Gets the actual section data\n",
    "\n",
    "      if kp_section.size < ir_padded.size: # Edge case for final section, as padded IR array may be longer than remaining data\n",
    "        np.pad(kp_section, (0, ir_padded.size - kp_section.size)) # Pads final section to length of padded IR\n",
    "\n",
    "      ir_scale_factor = self.df.at[i, \"kp_index\"] / 9\n",
    "      dry_scale_factor = 1 - ir_scale_factor # Creates a scale factor for dry signal, which is the inverse of ir_scale_factor\n",
    "      ir_padded_scaled = ir_padded * ir_scale_factor # Scaling the IR based on Kp index\n",
    "\n",
    "      Kp_Section = scipy.fft.fft(kp_section) # FFT of the kp_section\n",
    "      Ir_Padded_Scaled = scipy.fft.fft(ir_padded_scaled) # FFT of the padded scaled IR\n",
    "\n",
    "      Kp_Section_Reverb = Kp_Section * Ir_Padded_Scaled # Frequency domain multiplication\n",
    "\n",
    "      kp_section_reverb = np.real(scipy.fft.ifft(Kp_Section_Reverb)) # Returns the real section in the time domain\n",
    "\n",
    "      kp_section_mix = kp_section_reverb + (dry_scale_factor*kp_section) # Adding some of the dry signal back in for makeup gain\n",
    "\n",
    "      kp_section_mix_windowed = kp_section_mix * section_window # Windows the reverbed section\n",
    "\n",
    "      for sample_index, sample in enumerate(kp_section_mix_windowed): # Loops through mixed section\n",
    "        output_index = i * grain_len + sample_index # Corresponding index in output is the previous number of sections plus the current sample number\n",
    "        output_reverb[output_index] = sample # Adding mixed sample to output\n",
    "\n",
    "    output_reverb = self.__normalise(output_reverb)\n",
    "\n",
    "    rev_exec_time = time.time() - rev_start_time\n",
    "    print(\"Applied reverb to\", grains_type, \"grains in\", round(rev_exec_time, 3), \"secs\")\n",
    "\n",
    "    output_filepath = \"./output/\" + grains_type + \"_grains_output.wav\"\n",
    "    print(\"Writing output to:\", output_filepath)\n",
    "    sf.write(output_filepath, output_reverb, self.sr)\n",
    "    \n",
    "    return output_reverb\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  def reverb_both(self, path):\n",
    "    \"\"\"Utility method to apply reverb to both audios.\"\"\"\n",
    "    self.dp_synth_reverb = self.__apply_reverb(path, \"dp\")\n",
    "    self.midi_synth_reverb = self.__apply_reverb(path, \"midi\")\n",
    "\n",
    "    return (self.dp_synth_reverb, self.midi_synth_reverb)\n",
    "\n",
    "\n",
    "  def scatterplot(self, filetype='dp'):\n",
    "        \n",
    "    try:\n",
    "        if filetype == 'dp':\n",
    "            columns = ['dp_spec_band_grain_id', 'dp_rms_grain_id', 'dp_mfcc_grain_id']\n",
    "            title = 'Datapoint-segmented Synthesis'\n",
    "            grain_len = self.dp_grain_len\n",
    "            filename = './figures/dp_rolling_average_grain_id_distribution'\n",
    "        \n",
    "        elif filetype == 'midi':\n",
    "            columns = ['midi_spec_band_grain_id', 'midi_rms_grain_id', 'midi_mfcc_grain_id']\n",
    "            title = 'Midi-segmented Synthesis'\n",
    "            grain_len = self.midi_grain_len\n",
    "            filename = './figures/midi_rolling_average_grain_id_distribution'\n",
    "            \n",
    "        else: \n",
    "            raise ValueError() \n",
    "            \n",
    "    except ValueError:\n",
    "        print('invalid filetype')\n",
    "        \n",
    "        return \n",
    "    \n",
    "    plt.figure(figsize=(12,3))\n",
    "    plt.scatter(x=self.df.index*grain_len/self.sr, y=self.df[columns[0]], c='r', s=5, alpha=0.3, label=\"Spectral bandwidth grains\")\n",
    "    plt.scatter(x=self.df.index*grain_len/self.sr, y=self.df[columns[1]], c='g', s=5, alpha=0.3, label=\"RMS grains\")\n",
    "    plt.scatter(x=self.df.index*grain_len/self.sr, y=self.df[columns[2]], c='b', s=5, alpha=0.3, label=\"MFCC grains\")\n",
    "    plt.plot(self.df.index*grain_len/self.sr, self.df[columns].mean(axis=1),\n",
    "             c='black', alpha=0.7, lw=0.8, label=\"Average Grain ID Distribution\")\n",
    "    \n",
    "    plt.plot(self.df.index[::12]*grain_len/self.sr, \n",
    "             self.df[columns].mean(axis=1).rolling(12,step=12).mean(), c='white', label=\"Rolling Grain ID Average\")\n",
    "    plt.xlabel('Time (s)')\n",
    "    plt.ylabel('Grain ID')\n",
    "    plt.title(title)\n",
    "    plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))\n",
    "    plt.margins(0,0)\n",
    "    plt.tight_layout()\n",
    "    plt.savefig(filename, dpi=800)\n",
    "    plt.show()\n",
    "    \n",
    "   \n",
    "    \n",
    "  def __spectrogram(self, filetype='dp'):\n",
    "    \n",
    "    try:\n",
    "        if filetype == 'dp':\n",
    "            signal = self.dp_synth_reverb\n",
    "            title = 'Datapoint-segmented Synthesis'\n",
    "        \n",
    "        elif filetype == 'midi':\n",
    "            signal = self.midi_synth_reverb\n",
    "            title = 'Midi-segmented Synthesis'\n",
    "            \n",
    "        else: \n",
    "            raise ValueError() \n",
    "            \n",
    "    except ValueError:\n",
    "        print('invalid filetype')\n",
    "        \n",
    "        return \n",
    "    \n",
    "    \n",
    "    fig, ax = plt.subplots(figsize=(12,3))\n",
    "    data = librosa.display.specshow(librosa.amplitude_to_db(np.abs(librosa.stft(signal)), ref=np.max), \n",
    "                                    sr=self.sr, y_axis='log', x_axis='time', ax=ax)\n",
    "    \n",
    "    ax.set_title(title)\n",
    "    ax.set_ylabel('Frequency [Hz]')\n",
    "    ax.set_xlabel('Time [min]')\n",
    "    fig.colorbar(data, ax=ax, format=\"%+2.0f dB\") #format applies steps of 10dB to the colorbar \n",
    "    plt.show()\n",
    "    \n",
    "  def spectrograms(self):\n",
    "        \n",
    "    self.__spectrogram(filetype='dp')\n",
    "    self.__spectrogram(filetype='midi')\n",
    "        \n",
    "    print('--------Compiling video. Second coffee break maybe?--------')\n",
    "        \n",
    "        \n",
    "    \n",
    "    \n",
    "  def animate_graph(self):\n",
    "    \"\"\"Animation function for our dataset over 2016 datapoints. Saves the video in your directory once the cell finished running.\"\"\"\n",
    "    \n",
    "    self.df_anim = self.df.copy()\n",
    "    self.df_anim = self.df_anim[[\"proton_density\", \"ion_temp\", \"speed\", \"bz\", \"kp_index\"]]\n",
    "   \n",
    "    # Scaling data to a range of -1 to 1 for presentation\n",
    "    colnames = list(self.df_anim.columns)\n",
    "\n",
    "    # Scaling values\n",
    "    for column in colnames:\n",
    "        self.df_anim[column] = self.df_anim[column] / np.max(self.df_anim[column])\n",
    "\n",
    "    # Smoothing kp-index graph with interpolate function\n",
    "    self.kp_freq = 36\n",
    "    step_index = np.arange(0, self.dataset_len, self.kp_freq)\n",
    "    y = self.df_anim.iloc[step_index]['kp_index']\n",
    "    f = scipy.interpolate.make_interp_spline(step_index, y, k=2)\n",
    "    kp_index_smoothed = f(self.df.index)\n",
    "    self.df_anim['kp_index'] = kp_index_smoothed\n",
    "    \n",
    "    # Further scaling of Bz kp_index and speed values \n",
    "    self.df_anim[\"bz\"] = self.df_anim[\"bz\"]*0.5\n",
    "    self.df_anim['kp_index'] = self.df_anim['kp_index'] -0.5\n",
    "    self.df_anim[\"speed\"] = self.df_anim[\"speed\"] - self.df_anim[\"speed\"].min()\n",
    "    \n",
    "    # Creating empty lists for writing in values by animate function\n",
    "    x = []\n",
    "    proton_density = []\n",
    "    ion_temp = []\n",
    "    speed = []\n",
    "    bz = []\n",
    "    kp_index = []\n",
    "\n",
    "    # Counter variable for indexing\n",
    "    iterator = itertools.takewhile(lambda frame: frame < 2016, itertools.count(1))\n",
    "\n",
    "    # Figure and axes\n",
    "    fig = plt.figure(figsize=(16, 12), facecolor='black')\n",
    "    ax = plt.subplot(frameon=False) \n",
    "    \n",
    "    # Define ticks\n",
    "    means = self.df_anim.mean(axis=0).to_numpy()\n",
    "    offsets = np.arange(0, 1.5, 0.3)\n",
    "    ticks = means + offsets\n",
    "\n",
    "    def animate(i):\n",
    "        # Animation function which is getting called by FuncAnimation\n",
    "        \n",
    "        \n",
    "        # Writing from the dataset to lists\n",
    "        try:\n",
    "            next_iter = next(iterator)\n",
    "            x.append(next_iter)\n",
    "    \n",
    "            proton_density.append(self.df_anim['proton_density'].iat[next_iter] + offsets[0])\n",
    "            ion_temp.append((self.df_anim['ion_temp'].iat[next_iter]) + offsets[1])\n",
    "            speed.append((self.df_anim['speed'].iat[next_iter]) + offsets[2])\n",
    "            bz.append((self.df_anim['bz'].iat[next_iter]) + offsets[3])\n",
    "            kp_index.append((self.df_anim['kp_index'].iat[next_iter]) + offsets[4])\n",
    "    \n",
    "            #plotting data \n",
    "            plt.cla()\n",
    "    \n",
    "            ax.set_yticks(ticks,['Proton Density','Ion Temperature','Speed','BZ', 'KP Index'])\n",
    "            ax.get_yticklabels()[0].set_color('#00CF52')\n",
    "            ax.get_yticklabels()[1].set_color('#FF33EC')\n",
    "            ax.get_yticklabels()[2].set_color('#B53DFF')\n",
    "            ax.get_yticklabels()[3].set_color('#00C690')\n",
    "            ax.get_yticklabels()[4].set_color('#88F0F1')\n",
    "            ax.text(0.45, 1.0, 'AURORA',transform=ax.transAxes,ha=\"right\",\n",
    "                    color= '#00C690', fontweight='bold', fontsize=42)\n",
    "            ax.text(0.65, 1.0, \"VISUALIS\",transform=ax.transAxes,ha=\"right\",\n",
    "                    color='#88F0F1', fontweight='light',fontsize=42)\n",
    "           \n",
    "            \n",
    "            plt.plot(x,proton_density, c='#00CF52' )\n",
    "            plt.plot(x,ion_temp, c='#FF33EC')\n",
    "            plt.plot(x,speed, c='#B53DFF')\n",
    "            plt.plot(x,bz, c='#00C690')\n",
    "            plt.plot(x,kp_index, c='#88F0F1')\n",
    "    \n",
    "            plt.ylim(0,2.5)\n",
    "            \n",
    "            plt.tight_layout()\n",
    "        \n",
    "       \n",
    "            \n",
    "        # Catches the final iteration of the loop within FuncAnim\n",
    "        except StopIteration:\n",
    "            pass\n",
    "  \n",
    "    animation = anim.FuncAnimation(plt.gcf(), animate, interval=((self.dp_grain_len/self.sr) * 1000),\n",
    "                                           save_count=self.dataset_len)\n",
    "    \n",
    "    # Saving the animated graph and combining video and audio with moviepy\n",
    "    animation.save(\"./output/data_animated.gif\")\n",
    "    audioclip = AudioFileClip(\"./output/dp_grains_output.wav\")\n",
    "    videoclip = VideoFileClip(\"./output/data_animated.gif\")\n",
    "    videoclip_final = videoclip.set_audio(audioclip)\n",
    "    \n",
    "    # We only run the data visualisation for the datapoint-segmented synthesis - \n",
    "    # otherwise the program would run for 30 minutes\n",
    "    \n",
    "    videoclip_final.write_videofile(\"./output/dp_animation_synthesized.mp4\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Run this cell to compile output\n",
    "\n",
    "sonify = Sonify()\n",
    "\n",
    "# Read dataset from csv file\n",
    "sonify.read_data(\"./data/solar_wind_data_2003-10-27 - 2003-11-02 1min.xlsx\", \"./data/cleaned_data.csv\")\n",
    "\n",
    "# Read audio from corpus\n",
    "sonify.read_audio(\"./corpus/02_Dido White Flag.wav\")\n",
    "\n",
    "# Read MIDI from corpus\n",
    "sonify.read_midi(\"./corpus/02_Dido White Flag_adjusted_2.mid\")\n",
    "\n",
    "# Granulate the audio by both datapoint and MIDI segmentation, and create dataframes of grain IDs and start points\n",
    "sonify.build_grains_dataframes()\n",
    "\n",
    " # Map features from the grain dataframes to the dataset by copying across grain IDs\n",
    "sonify.map_features()\n",
    "\n",
    " # Scale the features used for processing parameters to suitable ranges\n",
    "sonify.scale_features()\n",
    "\n",
    "# Compound grains and apply mapping\n",
    "synth_outputs = sonify.synthesise_both()\n",
    "\n",
    "# # Apply the reverb\n",
    "outputs_reverb = sonify.reverb_both(\"./ir/ir_sydney_cathedral.wav\")\n",
    "\n",
    "# # Display the outputs of both segmentation methods\n",
    "display(Audio(outputs_reverb[0], rate=sonify.sr))\n",
    "display(Audio(outputs_reverb[1], rate=sonify.sr))\n",
    "\n",
    "# Visualisations\n",
    "sonify.scatterplot('dp')\n",
    "sonify.scatterplot('midi')\n",
    "sonify.spectrograms()\n",
    "\n",
    "# Aurora Visualis\n",
    "sonify.animate_graph()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.12"
  },
  "vscode": {
   "interpreter": {
    "hash": "33e85536a753e4289d1e3ca1904fea47091c61e030b0b2f85895695a2354c09d"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}