{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.7.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":19596,"databundleVersionId":1292430,"sourceType":"competition"}],"dockerImageVersionId":29962,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"![](https://64.media.tumblr.com/39e7a0085e9c2a5e87f70f89e259147c/tumblr_ow23gvvY2F1wrt8hqo1_1280.jpg)\nCan't help it..","metadata":{}},{"cell_type":"markdown","source":"# 🦜Cornell BirdSong Recognition","metadata":{}},{"cell_type":"markdown","source":"These are kind of my Notes from series:\n- [Audio Signal Processing](http://https://www.youtube.com/watch?v=bnHHVo3j124) \n- [MIR Website](http://musicinformationretrieval.com)","metadata":{}},{"cell_type":"markdown","source":"## Table of contents 📋\n1) Audio Signal Processing\n   - Mechanical Wave\n   - Information in Waves\n     - Frequency\n     - Pitch \n     - Timbre\n     - Sound Power,Sound Intensity\n     - Thresholds(of pain,Hearing)\n     - Loudness\n     - Intensity Level\n     - Complex Sound\n2) Librosa","metadata":{}},{"cell_type":"markdown","source":"## Audio Signal Processing\nProblem: How to use Audio for the Audio Classification??","metadata":{}},{"cell_type":"markdown","source":"What is a Audio or Sound??\n- Produced by vibration of an Object\n- Vibrations cause air molecules to oscilatte\n- change in air pressure creates a wave\n\nWhat is a mechanical wave(sound is a mechanical wave)?\n- Oscillation that travels through space\n- Energy travels from one point to another\n- the medium is deformed\n![](https://cdn.britannica.com/20/4320-050-9D00EC56/representations-Air-sound-wave-equilibrium-compressions-absence.jpg)","metadata":{}},{"cell_type":"markdown","source":"What information a waveform contains??\n- Frequency\n- Intesity\n- Timbre\nTypes of Waveforms:\n![](https://i.imgur.com/ubtrD6O.png)\nLarger ampltitude-> Louder\nLarger Period-> small frequency\n\nHearing Range of various animals:\n![](https://upload.wikimedia.org/wikipedia/commons/5/5d/Animal_hearing_frequency_range.svg)\n![](https://i.imgur.com/fvHJFyp.png)","metadata":{}},{"cell_type":"markdown","source":"Pitch:\n- We see frequency not in linear way but rather logarthmic perception\n- 2 frequencies are percieved to be almost similiarly if they differ by a power of 2\nMapping pitch to frequency:\n\n![](https://media.springernature.com/full/springer-static/image/art%3A10.1038%2Fs41598-017-18150-y/MediaObjects/41598_2017_18150_Fig1_HTML.jpg?as=webp)","metadata":{}},{"cell_type":"markdown","source":"Notes and Frequency\n![](https://newt.phys.unsw.edu.au/jw/graphics/notes.GIF)","metadata":{}},{"cell_type":"markdown","source":"Sound Power\n- Rate at which energy is transferred\n-  Energy per unit of time emitted by a sound source in all directions\n- Measured in watt(W)\n\nSound Intensity\n- Sound Power per unit area\n- Measured in W/m^2\n\nThresholds\n- human can percieve sounds with very small intensities\n- Threshold of hearing= 10^-12 W/m^-2\n- Threshold of Pain= 10 W/m^-2\n\nIntensity Level\n- Logaritmic Scale\n- Measured in decibels(dB)\n- Ratio between intensity Values\n- use an intesity of reference(TOH)\n- Every ~3Dbs\n![](http://hyperphysics.phy-astr.gsu.edu/hbase/Sound/imgsou/intens4.gif)\nLoudness\n- Subjective perception of sound intensity\n- depends on duration/frquency of a sound\n- Depends on age\n- Measured in phons\n![](https://upload.wikimedia.org/wikipedia/commons/thumb/4/47/Lindos1.svg/1024px-Lindos1.svg.png)","metadata":{}},{"cell_type":"markdown","source":"Timbre\n- Kind of color of sound\n- Difference between 2 seconds with same intensity,frequency,duration\n- Described with words like: bright,dark,dull,harsh,warm\nFeatures of timbre:\n-It is multidimensional\n- sound envelope\n![](https://www.audiolabs-erlangen.de/resources/MIR/FMP/data/C1/FMP_C1_F22a-23.png)\n- Harmonic content\n- Ampltitude/frquency mode","metadata":{}},{"cell_type":"code","source":"#Sound Envelope\n#- Attack-decay-sustain-release-model\nfrom IPython.display import IFrame\nIFrame('https://tonejs.github.io/examples/envelope.html', width=700, height=350)","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:58:53.177669Z","iopub.execute_input":"2023-07-17T15:58:53.178082Z","iopub.status.idle":"2023-07-17T15:58:53.187047Z","shell.execute_reply.started":"2023-07-17T15:58:53.178034Z","shell.execute_reply":"2023-07-17T15:58:53.18593Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Complex Sound\n- Superposition of sound\n- a partial is a sinusoid used to describe a sound\n- lowest partial is called fudamental frquency\n- harmonic partial is a frquency that's multiple of fundamental frquency\n- inharmonicity  indicates a deviation from harmonic partial\n- pitched instruments tend to be harmonic and percussive instruments tend to be inharmonic","metadata":{}},{"cell_type":"markdown","source":"# Librosa\nLibrosa is a Python package for music and audio processing by Brian McFee. ","metadata":{}},{"cell_type":"code","source":"pwd()","metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","execution":{"iopub.status.busy":"2023-07-17T15:58:53.300019Z","iopub.execute_input":"2023-07-17T15:58:53.300384Z","iopub.status.idle":"2023-07-17T15:58:53.306623Z","shell.execute_reply.started":"2023-07-17T15:58:53.300354Z","shell.execute_reply":"2023-07-17T15:58:53.305528Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cd ../input/birdsong-recognition/train_audio/aldfly","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:58:53.340129Z","iopub.execute_input":"2023-07-17T15:58:53.340494Z","iopub.status.idle":"2023-07-17T15:58:53.354967Z","shell.execute_reply.started":"2023-07-17T15:58:53.340464Z","shell.execute_reply":"2023-07-17T15:58:53.353917Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ls","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:58:53.399622Z","iopub.execute_input":"2023-07-17T15:58:53.400062Z","iopub.status.idle":"2023-07-17T15:58:54.481537Z","shell.execute_reply.started":"2023-07-17T15:58:53.400024Z","shell.execute_reply":"2023-07-17T15:58:54.480368Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Reading Audio\nUse librosa.load to load an audio file into an audio array. Return both the audio array as well as the sample rate:","metadata":{}},{"cell_type":"code","source":"!pip install ffmpeg","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:58:54.484058Z","iopub.execute_input":"2023-07-17T15:58:54.484359Z","iopub.status.idle":"2023-07-17T15:59:04.65229Z","shell.execute_reply.started":"2023-07-17T15:58:54.484326Z","shell.execute_reply":"2023-07-17T15:59:04.65144Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#loading an audio file into an audio array\nimport librosa\nx, sr = librosa.load('XC134874.mp3')","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:59:04.653601Z","iopub.execute_input":"2023-07-17T15:59:04.653936Z","iopub.status.idle":"2023-07-17T15:59:09.00305Z","shell.execute_reply.started":"2023-07-17T15:59:04.6539Z","shell.execute_reply":"2023-07-17T15:59:09.001912Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#To display the length of the audio array and samplerate\nprint(x.shape)\nprint(sr)","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:59:09.004458Z","iopub.execute_input":"2023-07-17T15:59:09.004828Z","iopub.status.idle":"2023-07-17T15:59:09.011156Z","shell.execute_reply.started":"2023-07-17T15:59:09.004781Z","shell.execute_reply":"2023-07-17T15:59:09.009824Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Visualizing Audio","metadata":{}},{"cell_type":"code","source":"%matplotlib inline\nimport matplotlib.pyplot as plt\nimport librosa.display","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:59:09.015004Z","iopub.execute_input":"2023-07-17T15:59:09.015438Z","iopub.status.idle":"2023-07-17T15:59:09.025726Z","shell.execute_reply.started":"2023-07-17T15:59:09.015391Z","shell.execute_reply":"2023-07-17T15:59:09.024865Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#plotting the audio array with librosa.display.wavplot\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr=sr)","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:59:09.029021Z","iopub.execute_input":"2023-07-17T15:59:09.02944Z","iopub.status.idle":"2023-07-17T15:59:09.266029Z","shell.execute_reply.started":"2023-07-17T15:59:09.029396Z","shell.execute_reply":"2023-07-17T15:59:09.264847Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#displaying the spectrogram using librosa.display.specshow\nX = librosa.stft(x)\nXdb = librosa.amplitude_to_db(abs(X))\nplt.figure(figsize=(14, 5))\nlibrosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='hz')\nimport IPython.display as ipd\nipd.Audio('XC134874.mp3') # load a local WAV file","metadata":{"execution":{"iopub.status.busy":"2023-07-17T15:59:09.267876Z","iopub.execute_input":"2023-07-17T15:59:09.26827Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"we must extract the characteristics of our audio signal that are most relevant to the problem we are trying to solve.process is known as feature extraction.","metadata":{}},{"cell_type":"code","source":"%matplotlib inline\nfrom pathlib import Path\nimport numpy, scipy, matplotlib.pyplot as plt, sklearn, urllib, IPython.display as ipd\nimport librosa, librosa.display\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import warnings\nwarnings.filterwarnings(\"ignore\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"signals = [\n    librosa.load(p)[0] for p in Path().glob('XC*.mp3')\n]","metadata":{"execution":{"iopub.status.idle":"2023-07-17T16:04:08.280361Z","shell.execute_reply.started":"2023-07-17T15:59:10.214677Z","shell.execute_reply":"2023-07-17T16:04:08.278858Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"len(signals)","metadata":{"execution":{"iopub.status.busy":"2023-07-17T16:04:08.282253Z","iopub.execute_input":"2023-07-17T16:04:08.282605Z","iopub.status.idle":"2023-07-17T16:04:08.288527Z","shell.execute_reply.started":"2023-07-17T16:04:08.282566Z","shell.execute_reply":"2023-07-17T16:04:08.287801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(100, 500))\nfor i, x in enumerate(signals):\n    plt.subplot(25, 4, i+1)\n    librosa.display.waveplot(x[:10000])\n    plt.ylim(-1, 1)","metadata":{"execution":{"iopub.status.busy":"2023-07-17T16:04:08.289977Z","iopub.execute_input":"2023-07-17T16:04:08.29031Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#let's now construct a feature vector >feature vector is simply a collection of features.\ndef extract_features(signal):\n    return [\n        librosa.feature.zero_crossing_rate(signal)[0, 0],\n        librosa.feature.spectral_centroid(signal)[0, 0],\n    ]","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"features = numpy.array([extract_features(x) for x in signals])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nplt.hist(features[:,0], color='b', range=(0, 0.2), alpha=0.5, bins=20)\nplt.legend(('Aldfly signals'))\nplt.xlabel('Zero Crossing Rate')\nplt.ylabel('Count')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nplt.hist(features[:,1], color='b', range=(0, 4000), bins=30, alpha=0.6)\nplt.legend(('Aldfly Signals'))\nplt.xlabel('Spectral Centroid (frequency bin)')\nplt.ylabel('Count')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#this part does the feature scaling so that those fetures can worked on with each other\n#so like what we do here is bringing them all in range of -1,1 as neccessary","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"feature_table = numpy.vstack((features))\nprint(feature_table.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"scaler = sklearn.preprocessing.MinMaxScaler(feature_range=(-1, 1))\ntraining_features = scaler.fit_transform(feature_table)\nprint(training_features.min(axis=0))\nprint(training_features.max(axis=0))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.scatter(training_features[:50,0], training_features[:50,1], c='b')\nplt.scatter(training_features[50:,0], training_features[50:,1], c='r')\nplt.xlabel('Zero Crossing Rate')\nplt.ylabel('Spectral Centroid')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Segmentation\nIn audio processing, it is common to operate on one frame at a time using a constant frame size and hop size (i.e. increment). Frames are typically chosen to be 10 to 100 ms in duration.\n\nLet's create an audio signal consisting of a pure tone that gradually gets louder. Then, we will segment the signal and compute the root mean square (RMS) energy for each frame.\n\nFirst, set our parameters:","metadata":{}},{"cell_type":"code","source":"#Segmentation \n%matplotlib inline\nimport seaborn\nimport numpy, scipy, matplotlib.pyplot as plt, IPython.display as ipd\nimport librosa, librosa.display\nplt.rcParams['figure.figsize'] = (11, 5)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"T = 3.0      # duration in seconds\nsr = 22050   # sampling rate in Hertz\namplitude = numpy.logspace(-3, 0, int(T*sr), endpoint=False, base=10.0) # time-varying amplitude\nprint amplitude.min(), amplitude.max()\n# starts at 110 Hz, ends at 880 Hz\n#create the signal\nt = numpy.linspace(0, T, int(T*sr), endpoint=False)\nx = amplitude*numpy.sin(2*numpy.pi*440*t)\n\nipd.Audio(x, rate=sr)\n#Plot the signal:\nlibrosa.display.waveplot(x, sr=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Segmentation Using Python List Comprehensions\n#In Python, you can use a standard list comprehension to perform segmentation of a signal and compute RMSE at the same time.\n#Initialize segmentation parameters:\n\nframe_length = 1024\nhop_length = 512\n#Define a helper function:\n\ndef rmse(x):\n    return numpy.sqrt(numpy.mean(x**2))\n\n#Using a list comprehension, plot the RMSE for each frame on a log-y axis:\n\nplt.semilogy([rmse(x[i:i+frame_length])\n              for i in range(0, len(x), hop_length)])\n\n#librosa.util.frame\n#Given a signal, librosa.util.frame will produce a list of uniformly sized frames:\n\nframes = librosa.util.frame(x, frame_length=frame_length, hop_length=hop_length)\nplt.semilogy([rmse(frame) for frame in frames.T])\n\n#That being said, in librosa, manual segmentation of a signal is often unnecessary, because the feature extraction methods themselves do segmentation for you.","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Zero Crossing Rate\n#The zero crossing rate indicates the number of times that a signal crosses the horizontal axis.\n#Let's load a signal:\n\nx, sr = librosa.load('XC134874.mp3')\n#Listen to the signal:\nipd.Audio(x, rate=sr)\n\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Let's zoom in:\n\nn0 = 6500\nn1 = 7500\nplt.figure(figsize=(14, 5))\nplt.plot(x[n0:n1])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#I count five zero crossings. Let's compute the zero crossings using librosa.\n\nzero_crossings = librosa.zero_crossings(x[n0:n1], pad=False)\nzero_crossings.shape\n#That computed a binary mask where True indicates the presence of a zero crossing. To find the total number of zero crossings, use sum:\nprint(sum(zero_crossings))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#To find the zero-crossing rate over time, use zero_crossing_rate:\nzcrs = librosa.feature.zero_crossing_rate(x)\nprint(zcrs.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nplt.plot(zcrs[0])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Note how the high zero-crossing rate corresponds to the presence of the snare drum.\n\n#The reason for the high rate near the beginning is because the silence oscillates quietly around zero:\n\nplt.figure(figsize=(14, 5))\nplt.plot(x[:1000])\nplt.ylim(-0.0001, 0.0001)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#A simple hack around this is to add a small constant before computing the zero crossing rate:\n\nzcrs = librosa.feature.zero_crossing_rate(x + 0.0001)\nplt.figure(figsize=(14, 5))\nplt.plot(zcrs[0])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# FourierTransform\n","metadata":{}},{"cell_type":"code","source":"import seaborn\nimport numpy, scipy, matplotlib.pyplot as plt, librosa, IPython.display as ipd","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import urllib\nfilename = 'XC134874.mp3'\nx, sr = librosa.load(filename)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print x.shape\nprint sr","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Fourier Transform\nThe fourier transform is one of the most fundamental operations in applied mathematics and signal processing.\n\nIt transforms our time-domain signal into the frequency domain. Whereas the time domain expresses our signal as a sequence of samples, the frequency domain expresses our signal as a superposition of sinusoids of varying magnitudes, frequencies, and phase offsets.\n\nTo compute a Fourier transform in NumPy or SciPy, use scipy.fft:","metadata":{}},{"cell_type":"code","source":"X = scipy.fft(x)\nX_mag = numpy.absolute(X)\nf = numpy.linspace(0, sr, len(X_mag)) # frequency variable","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the spectrum:\n\nplt.figure(figsize=(13, 5))\nplt.plot(f, X_mag) # magnitude spectrum\nplt.xlabel('Frequency (Hz)')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Zoom in:\nplt.figure(figsize=(13, 5))\nplt.plot(f[:5000], X_mag[:5000])\nplt.xlabel('Frequency (Hz)')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Short Time Fourier Transform\n","metadata":{}},{"cell_type":"code","source":"#loadingfile\nx, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#librosa.stft computes a STFT. \n#We provide it a frame size, i.e. the size of the FFT, and a hop length, i.e. the frame increment:\nhop_length = 512\nn_fft = 2048\nX = librosa.stft(x, n_fft=n_fft, hop_length=hop_length)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#To convert the hop length and frame size to units of seconds:\n\nfloat(hop_length)/sr # units of seconds","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"float(n_fft)/sr  # units of seconds","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#For real-valued signals, the Fourier transform is symmetric about the midpoint. Therefore, librosa.stft only retains one half of the output:\nX.shape\n#This STFT has 1025 frequency bins and 9813 frames in time.","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectrogram\nwe often care about the spectral magnitude and not the phrase content\nThe spectrogram (Wikipedia; FMP, p. 29, 55) shows the the intensity of frequencies over time. A spectrogram is simply the squared magnitude of the STFT:\n\n$$ S(m, \\omega) = \\left| X(m, \\omega) \\right|^2 $$\n\nThe human perception of sound intensity is logarithmic in nature. Therefore, we are often interested in the log amplitude:","metadata":{}},{"cell_type":"code","source":"S = librosa.amplitude_to_db(abs(X))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To display any type of spectrogram in librosa, use librosa.display.specshow.","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(15, 5))\nlibrosa.display.specshow(S, sr=sr, hop_length=hop_length, x_axis='time', y_axis='linear')\nplt.colorbar(format='%+2.0f dB')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Mel-spectrogram\nlibrosa has some outstanding spectral representations, including librosa.feature.melspectrogram:","metadata":{}},{"cell_type":"code","source":"hop_length = 256\nS = librosa.feature.melspectrogram(x, sr=sr, n_fft=4096, hop_length=hop_length)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The human perception of sound intensity is logarithmic in nature. Therefore, like the STFT-based spectrogram, we are often interested in the log amplitude:\nTo display any type of spectrogram in librosa, use librosa.display.specshow.","metadata":{}},{"cell_type":"code","source":"logS = librosa.power_to_db(abs(S))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(15, 5))\nlibrosa.display.specshow(logS, sr=sr, hop_length=hop_length, x_axis='time', y_axis='mel')\nplt.colorbar(format='%+2.0f dB')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Using y_axis=mel plots the y-axis on the mel scale which is similar to the $\\log (1 + f)$ function:\n\n$$ m = 2595 \\log_{10} \\left(1 + \\frac{f}{700} \\right) $$","metadata":{}},{"cell_type":"markdown","source":"# librosa.cqt\nUnlike the Fourier transform, but similar to the mel scale, the constant-Q transform uses a logarithmically spaced frequency axis.\n\nTo plot a constant-Q spectrogram, will use librosa.cqt:","metadata":{}},{"cell_type":"code","source":"fmin = librosa.midi_to_hz(36)\nC = librosa.cqt(x, sr=sr, fmin=fmin, n_bins=72)\nlogC = librosa.amplitude_to_db(abs(C))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(15, 5))\nlibrosa.display.specshow(logC, sr=sr, x_axis='time', y_axis='cqt_note', fmin=fmin, cmap='coolwarm')\nplt.colorbar(format='%+2.0f dB')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Constant-Q Transform and Chroma\nConstant-Q Transform\nUnlike the Fourier transform, but similar to the mel scale, the constant-Q transform (Wikipedia) uses a logarithmically spaced frequency axis. For more information, read the original paper:\n\nJudith C. Brown, \"Calculation of a constant Q spectral transform,\" J. Acoust. Soc. Am., 89(1):425–434, 1991.\nLet's load a file:","metadata":{}},{"cell_type":"code","source":"ipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To compute a constant-Q spectrogram, will use librosa.cqt:","metadata":{}},{"cell_type":"code","source":"fmin = librosa.midi_to_hz(36)\nhop_length = 512\nC = librosa.cqt(x, sr=sr, fmin=fmin, n_bins=72, hop_length=hop_length)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Display:\nlogC = librosa.amplitude_to_db(numpy.abs(C))\nplt.figure(figsize=(15, 5))\nlibrosa.display.specshow(logC, sr=sr, x_axis='time', y_axis='cqt_note', fmin=fmin, cmap='coolwarm')\n#Note how each frequency bin corresponds to one MIDI pitch number.","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Chroma\nA chroma vector (Wikipedia) (FMP, p. 123) is a typically a 12-element feature vector indicating how much energy of each pitch class, {C, C#, D, D#, E, ..., B}, is present in the signal.\n\nlibrosa.feature.chroma_stft","metadata":{}},{"cell_type":"code","source":"chromagram = librosa.feature.chroma_stft(x, sr=sr, hop_length=hop_length)\nplt.figure(figsize=(15, 5))\nlibrosa.display.specshow(chromagram, x_axis='time', y_axis='chroma', hop_length=hop_length, cmap='coolwarm')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# librosa.feature.chroma_cqt","metadata":{}},{"cell_type":"code","source":"chromagram = librosa.feature.chroma_cqt(x, sr=sr, hop_length=hop_length)\nplt.figure(figsize=(15, 5))\nlibrosa.display.specshow(chromagram, x_axis='time', y_axis='chroma', hop_length=hop_length, cmap='coolwarm')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Chroma energy normalized statistics (CENS) (FMP, p. 375). The main idea of CENS features is that taking statistics over large windows smooths local deviations in tempo, articulation, and musical ornaments such as trills and arpeggiated chords. CENS are best used for tasks such as audio matching and similarity.\n\nlibrosa.feature.chroma_cens","metadata":{}},{"cell_type":"code","source":"chromagram = librosa.feature.chroma_cens(x, sr=sr, hop_length=hop_length)\nplt.figure(figsize=(15, 5))\nlibrosa.display.specshow(chromagram, x_axis='time', y_axis='chroma', hop_length=hop_length, cmap='coolwarm')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Magnitude Scaling\nOften, the raw amplitude of a signal in the time- or frequency-domain is not as perceptually relevant to humans as the amplitude converted into other units, e.g. using a logarithmic scale.\n\nFor example, let's consider a pure tone whose amplitude grows louder linearly. Define the time variable:","metadata":{}},{"cell_type":"code","source":"T = 4.0      # duration in seconds\nsr = 22050   # sampling rate in Hertz\nt = numpy.linspace(0, T, int(T*sr), endpoint=False)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Create a signal whose amplitude grows linearly:\namplitude = numpy.linspace(0, 1, int(T*sr), endpoint=False) # time-varying amplitude\nx = amplitude*numpy.sin(2*numpy.pi*440*t)\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"librosa.display.waveplot(x, sr=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Now consider a signal whose amplitude grows exponentially, i.e. the logarithm of the amplitude is linear:\n\namplitude = numpy.logspace(-2, 0, int(T*sr), endpoint=False, base=10.0)\nx = amplitude*numpy.sin(2*numpy.pi*440*t)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"librosa.display.waveplot(x, sr=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Even though the amplitude grows exponentially, to us, the increase in loudness seems more gradual. This phenomenon is an example of the Weber-Fechner law (Wikipedia) which states that the relationship between a stimulus and human perception is logarithmic.","metadata":{}},{"cell_type":"markdown","source":"# Spectrogram Visualization: Linear Amplitude¶\nLet's plot a magnitude spectrogram where the colorbar is a linear function of the spectrogram values, i.e. just plot the raw values.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3', duration=25)\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X = librosa.stft(x)\nX.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Raw amplitude:\nXmag = abs(X)\nlibrosa.display.specshow(Xmag, sr=sr, x_axis='time', y_axis='log')\nplt.colorbar()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectrogram Visualization: Log Amplitude\nNow let's plot a magnitude spectrogram where the colorbar is a logarithmic function of the spectrogram values.\n\nDecibel (Wikipedia)\n\nlibrosa.amplitude_to_db:","metadata":{}},{"cell_type":"code","source":"Xdb = librosa.amplitude_to_db(Xmag)\nlibrosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='log')\nplt.colorbar()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"One common variant is the  log(1+λx)  function, sometimes known as logarithmic compression (FMP, p. 125). This function operates like  y=λx  when  λx  is small, but it operates like  y=logλx  when  λx  is large.","metadata":{}},{"cell_type":"code","source":"Xmag = numpy.log10(1+10*abs(X))\nlibrosa.display.specshow(Xmag, sr=sr, x_axis='time', y_axis='log', cmap=\"gray_r\")\nplt.colorbar()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectrogram Visualization: Perceptual Weighting\nlibrosa.perceptual_weighting:","metadata":{}},{"cell_type":"code","source":"freqs = librosa.core.fft_frequencies(sr=sr)\nXmag = librosa.perceptual_weighting(abs(X)**2, freqs)\nlibrosa.display.specshow(Xmag, sr=sr, x_axis='time', y_axis='log')\nplt.colorbar()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"","metadata":{}},{"cell_type":"code","source":"%matplotlib inline\nimport seaborn\nimport numpy, scipy, matplotlib.pyplot as plt, IPython.display as ipd, sklearn\nimport librosa, librosa.display\nplt.rcParams['figure.figsize'] = (14, 5)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectral Features\nFor classification, we're going to be using new features in our arsenal: spectral moments (centroid, bandwidth, skewness, kurtosis) and other spectral statistics.\n\nMoments is a term used in physics and statistics. There are raw moments and central moments.\n\nYou are probably already familiar with two examples of moments: mean and variance. The first raw moment is known as the mean. The second central moment is known as the variance.\n\nSpectral Centroid\nLoad an audio file:","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The spectral centroid (Wikipedia) indicates at which frequency the energy of a spectrum is centered upon. This is like a weighted mean:\n\nfc=∑kS(k)f(k)∑kS(k)\n \nwhere  S(k)  is the spectral magnitude at frequency bin  k ,  f(k)  is the frequency at bin  k .\n\nlibrosa.feature.spectral_centroid computes the spectral centroid for each frame in a signal:","metadata":{}},{"cell_type":"code","source":"spectral_centroids = librosa.feature.spectral_centroid(x, sr=sr)[0]\nspectral_centroids.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Compute the time variable for visualization:\n\nframes = range(len(spectral_centroids))\nt = librosa.frames_to_time(frames)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Define a helper function to normalize the spectral centroid for visualization:","metadata":{}},{"cell_type":"code","source":"def normalize(x, axis=0):\n    return sklearn.preprocessing.minmax_scale(x, axis=axis)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Plot the spectral centroid along with the waveform:","metadata":{}},{"cell_type":"code","source":"librosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_centroids), color='r') # normalize for visualization purposes","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\nSimilar to the zero crossing rate, there is a spurious rise in spectral centroid at the beginning of the signal. That is because the silence at the beginning has such small amplitude that high frequency components have a chance to dominate. One hack around this is to add a small constant before computing the spectral centroid, thus shifting the centroid toward zero at quiet portions:","metadata":{}},{"cell_type":"code","source":"spectral_centroids = librosa.feature.spectral_centroid(x+0.01, sr=sr)[0]\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_centroids), color='r') # normalize for visualization purposes","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectral Bandwidth\nlibrosa.feature.spectral_bandwidth computes the order- p  spectral bandwidth:\n\n(∑kS(k)(f(k)−fc)p)1p\n \nwhere  S(k)  is the spectral magnitude at frequency bin  k ,  f(k)  is the frequency at bin  k , and  fc  is the spectral centroid. When  p=2 , this is like a weighted standard deviation.","metadata":{}},{"cell_type":"code","source":"spectral_bandwidth_2 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr)[0]\nspectral_bandwidth_3 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr, p=3)[0]\nspectral_bandwidth_4 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr, p=4)[0]\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_bandwidth_2), color='r')\nplt.plot(t, normalize(spectral_bandwidth_3), color='g')\nplt.plot(t, normalize(spectral_bandwidth_4), color='y')\nplt.legend(('p = 2', 'p = 3', 'p = 4'))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectral Contrast\nSpectral contrast considers the spectral peak, the spectral valley, and their difference in each frequency subband. For more information:\n\nJiang, Dan-Ning, Lie Lu, Hong-Jiang Zhang, Jian-Hua Tao, and Lian-Hong Cai. “Music type classification by spectral contrast feature.” In Multimedia and Expo, 2002. ICME‘02. Proceedings. 2002 IEEE International Conference on, vol. 1, pp. 113-116. IEEE, 2002.\nlibrosa.feature.spectral_contrast computes the spectral contrast for six subbands for each time frame:","metadata":{}},{"cell_type":"code","source":"spectral_contrast = librosa.feature.spectral_contrast(x, sr=sr)\nspectral_contrast.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.imshow(normalize(spectral_contrast, axis=1), aspect='auto', origin='lower', cmap='coolwarm')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectral Rolloff\nSpectral rolloff is the frequency below which a specified percentage of the total spectral energy, e.g. 85%, lies.\n\nlibrosa.feature.spectral_rolloff computes the rolloff frequency for each frame in a signal:","metadata":{}},{"cell_type":"code","source":"spectral_rolloff = librosa.feature.spectral_rolloff(x+0.01, sr=sr)[0]\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_rolloff), color='r')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Autocorrelation\nThe autocorrelation of a signal describes the similarity of a signal against a time-shifted version of itself. For a signal $x$, the autocorrelation $r$ is:\n\n$$ r(k) = \\sum_n x(n) x(n-k) $$\n\nIn this equation, $k$ is often called the lag parameter. $r(k)$ is maximized at $k = 0$ and is symmetric about $k$.\n\nThe autocorrelation is useful for finding repeated patterns in a signal. For example, at short lags, the autocorrelation can tell us something about the signal's fundamental frequency. For longer lags, the autocorrelation may tell us something about the tempo of a musical signal.\n\nLet's load a file:","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# numpy.correlate\nThere are two ways we can compute the autocorrelation in Python. The first method is numpy.correlate:","metadata":{}},{"cell_type":"code","source":"# Because the autocorrelation produces a symmetric signal, we only care about the \"right half\".\nr = numpy.correlate(x, x, mode='full')[len(x)-1:]\nprint(x.shape, r.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the autocorrelation:\nplt.figure(figsize=(14, 5))\nplt.plot(r[:10000])\nplt.xlabel('Lag (samples)')\nplt.xlim(0, 10000)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# librosa.autocorrelate\nThe second method is librosa.autocorrelate:","metadata":{}},{"cell_type":"code","source":"r = librosa.autocorrelate(x, max_size=10000)\nprint(r.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nplt.plot(r)\nplt.xlabel('Lag (samples)')\nplt.xlim(0, 10000)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"librosa.autocorrelate conveniently only keeps one half of the autocorrelation function, since the autocorrelation is symmetric. Also, the max_size parameter prevents unnecessary calculations.\n\n# Pitch Estimation\nThe autocorrelation is used to find repeated patterns within a signal. For musical signals, a repeated pattern can correspond to a pitch period. We can therefore use the autocorrelation function to estimate the pitch in a musical signal.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Compute and plot the autocorrelation:\nr = librosa.autocorrelate(x, max_size=5000)\nplt.figure(figsize=(14, 5))\nplt.plot(r[:200])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The autocorrelation always has a maximum at zero, i.e. zero lag. We want to identify the maximum outside of the peak centered at zero. Therefore, we might choose only to search within a range of reasonable pitches:","metadata":{}},{"cell_type":"code","source":"midi_hi = 120.0\nmidi_lo = 12.0\nf_hi = librosa.midi_to_hz(midi_hi)\nf_lo = librosa.midi_to_hz(midi_lo)\nt_lo = sr/f_hi\nt_hi = sr/f_lo","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(f_lo, f_hi)\nprint(t_lo, t_hi)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Set invalid pitch candidates to zero:","metadata":{}},{"cell_type":"code","source":"r[:int(t_lo)] = 0\nr[int(t_hi):] = 0\nplt.figure(figsize=(14, 5))\nplt.plot(r[:1400])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Pitch Transcription Exercise\nLoad an audio file.","metadata":{}},{"cell_type":"code","source":"%matplotlib inline\nimport seaborn\nimport numpy, scipy, IPython.display as ipd, matplotlib.pyplot as plt\nimport librosa, librosa.display\nplt.rcParams['figure.figsize'] = (14, 5)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"filename = 'XC134874.mp3'\nx, sr = librosa.load(filename)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Play the audio file.\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Display the CQT of the signal.\nbins_per_octave = 36\ncqt = librosa.cqt(x, sr=sr, n_bins=300, bins_per_octave=bins_per_octave)\nlog_cqt = librosa.logamplitude(cqt)\ncqt.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"librosa.display.specshow(log_cqt, sr=sr, x_axis='time', y_axis='cqt_note', \n                         bins_per_octave=bins_per_octave)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Goal: Identify the pitch of each note and replace each note with a pure tone of that pitch.","metadata":{}},{"cell_type":"markdown","source":"# Step 1: Detect Onsets¶\nTo accurately detect onsets, it may be helpful to see what the novelty function looks like:","metadata":{}},{"cell_type":"code","source":"hop_length = 100\nonset_env = librosa.onset.onset_strength(x, sr=sr, hop_length=hop_length)\nplt.plot(onset_env)\nplt.xlim(0, len(onset_env))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Among the obvious large peaks, there are many smaller peaks. We want to choose parameters which preserve the large peaks while ignoring the small peaks.\n\nNext, we try to detect onsets. For more details, see librosa.onset.onset_detect and librosa.util.peak_pick.","metadata":{}},{"cell_type":"code","source":"onset_samples = librosa.onset.onset_detect(x,\n                                           sr=sr, units='samples', \n                                           hop_length=hop_length, \n                                           backtrack=False,\n                                           pre_max=20,\n                                           post_max=20,\n                                           pre_avg=100,\n                                           post_avg=100,\n                                           delta=0.2,\n                                           wait=0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"onset_samples","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Let's pad the onsets with the beginning and end of the signal.\nonset_boundaries = numpy.concatenate([[0], onset_samples, [len(x)]])\nprint(onset_boundaries)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Convert the onsets to units of seconds:\nonset_times = librosa.samples_to_time(onset_boundaries, sr=sr)\nonset_times","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Display the results of the onset detection:\nlibrosa.display.waveplot(x, sr=sr)\nplt.vlines(onset_times, -1, 1, color='r')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Step 2: Estimate Pitch\nEstimate pitch using the autocorrelation method:","metadata":{}},{"cell_type":"code","source":"def estimate_pitch(segment, sr, fmin=50.0, fmax=2000.0):\n    \n    # Compute autocorrelation of input segment.\n    r = librosa.autocorrelate(segment)\n    \n    # Define lower and upper limits for the autocorrelation argmax.\n    i_min = sr/fmax\n    i_max = sr/fmin\n    r[:int(i_min)] = 0\n    r[int(i_max):] = 0\n    \n    # Find the location of the maximum autocorrelation.\n    i = r.argmax()\n    f0 = float(sr)/i\n    return f0","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Step 3: Generate Pure Tone\nCreate a function to generate a pure tone at the specified frequency:","metadata":{}},{"cell_type":"code","source":"def generate_sine(f0, sr, n_duration):\n    n = numpy.arange(n_duration)\n    return 0.2*numpy.sin(2*numpy.pi*f0*n/float(sr))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Step 4: Put it together\nCreate a helper function for use in a list comprehension:","metadata":{}},{"cell_type":"code","source":"def estimate_pitch_and_generate_sine(x, onset_samples, i, sr):\n    n0 = onset_samples[i]\n    n1 = onset_samples[i+1]\n    f0 = estimate_pitch(x[n0:n1], sr)\n    return generate_sine(f0, sr, n1-n0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Use a list comprehension to concatenate the synthesized segments:","metadata":{}},{"cell_type":"code","source":"y = numpy.concatenate([\n    estimate_pitch_and_generate_sine(x, onset_boundaries, i, sr=sr)\n    for i in range(len(onset_boundaries)-1)\n])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Play the synthesized transcription.\nipd.Audio(y, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the CQT of the synthesized transcription.\ncqt = librosa.cqt(y, sr=sr)\nlibrosa.display.specshow(abs(cqt), sr=sr, x_axis='time', y_axis='cqt_note')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Novelty Functions\nTo detect note onsets, we want to locate sudden changes in the audio signal that mark the beginning of transient regions. Often, an increase in the signal's amplitude envelope will denote an onset candidate. However, that is not always the case, for notes can change from one pitch to another without changing amplitude, e.g. a violin playing slurred notes.\n\nNovelty functions are functions which denote local changes in signal properties such as energy or spectral content. We will look at two novelty functions:\n\nEnergy-based novelty functions (FMP, p. 306)\nSpectral-based novelty functions (FMP, p. 309)","metadata":{}},{"cell_type":"markdown","source":"Energy-based Novelty Functions\nPlaying a note often coincides with a sudden increase in signal energy. To detect this sudden increase, we will compute an energy novelty function (FMP, p. 307):\n\nCompute the short-time energy in the signal.\nCompute the first-order difference in the energy.\nHalf-wave rectify the first-order difference.\nFirst, load an audio file into the NumPy array x and sampling rate sr.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nprint(x.shape, sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the signal:\n\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr)\n#Listen:\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# RMS Energy\nlibrosa.feature.rmse returns the root-mean-square (RMS) energy for each frame of audio. We will compute the RMS energy as well as its first-order difference.","metadata":{}},{"cell_type":"code","source":"hop_length = 512\nframe_length = 1024\nrmse = librosa.feature.rmse(x, frame_length=frame_length, hop_length=hop_length).flatten()\nrmse_diff = numpy.zeros_like(rmse)\nrmse_diff[1:] = numpy.diff(rmse)\nprint(rmse.shape)\nprint(rmse_diff.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To obtain an energy novelty function, we perform half-wave rectification (FMP, p. 307) on rmse_diff, i.e. any negative values are set to zero. Equivalently, we can apply the function  max(0,x) :","metadata":{}},{"cell_type":"code","source":"energy_novelty = numpy.max([numpy.zeros_like(rmse_diff), rmse_diff], axis=0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot all three functions together:\n\nframes = numpy.arange(len(rmse))\nt = librosa.frames_to_time(frames, sr=sr)\nplt.figure(figsize=(15, 6))\nplt.plot(t, rmse, 'b--', t, rmse_diff, 'g--^', t, energy_novelty, 'r-')\nplt.xlim(0, t.max())\nplt.xlabel('Time (sec)')\nplt.legend(('RMSE', 'delta RMSE', 'energy novelty')) ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Log Energy\nThe human perception of sound intensity is logarithmic in nature. To account for this property, we can apply a logarithm function to the energy before taking the first-order difference.\n\nBecause  log(x)  diverges as  x  approaches zero, a common alternative is to use  log(1+λx) . This function equals zero when  x  is zero, but it behaves like  log(λx)  when  λx  is large. This operation is sometimes called logarithmic compression (FMP, p. 310).","metadata":{}},{"cell_type":"code","source":"log_rmse = numpy.log1p(10*rmse)\nlog_rmse_diff = numpy.zeros_like(log_rmse)\nlog_rmse_diff[1:] = numpy.diff(log_rmse)\nlog_energy_novelty = numpy.max([numpy.zeros_like(log_rmse_diff), log_rmse_diff], axis=0)\nplt.figure(figsize=(15, 6))\nplt.plot(t, log_rmse, 'b--', t, log_rmse_diff, 'g--^', t, log_energy_novelty, 'r-')\nplt.xlim(0, t.max())\nplt.xlabel('Time (sec)')\nplt.legend(('log RMSE', 'delta log RMSE', 'log energy novelty')) ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spectral-based Novelty Functions¶\nThere are two problems with the energy novelty function:\n\nIt is sensitive to energy fluctuations belonging to the same note.\nIt is not sensitive to spectral fluctuations between notes where amplitude remains the same.\nFor example, consider the following audio signal composed of pure tones of equal magnitude:","metadata":{}},{"cell_type":"code","source":"sr = 22050\ndef generate_tone(midi):\n    T = 0.5\n    t = numpy.linspace(0, T, int(T*sr), endpoint=False)\n    f = librosa.midi_to_hz(midi)\n    return numpy.sin(2*numpy.pi*f*t)\nx = numpy.concatenate([generate_tone(midi) for midi in [48, 52, 55, 60, 64, 67, 72, 76, 79, 84]])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Listen:\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#The energy novelty function remains roughly constant:\nhop_length = 512\nframe_length = 1024\nrmse = librosa.feature.rmse(x, frame_length=frame_length, hop_length=hop_length).flatten()\nrmse_diff = numpy.zeros_like(rmse)\nrmse_diff[1:] = numpy.diff(rmse)\nenergy_novelty = numpy.max([numpy.zeros_like(rmse_diff), rmse_diff], axis=0)\nframes = numpy.arange(len(rmse))\nt = librosa.frames_to_time(frames, sr=sr)\nplt.figure(figsize=(15, 4))\nplt.plot(t, rmse, 'b--', t, rmse_diff, 'g--^', t, energy_novelty, 'r-')\nplt.xlim(0, t.max())\nplt.xlabel('Time (sec)')\nplt.legend(('RMSE', 'delta RMSE', 'energy novelty')) ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Instead, we will compute a spectral novelty function (FMP, p. 309):\n\nCompute the log-amplitude spectrogram.\nWithin each frequency bin,  k , compute the energy novelty function as shown earlier, i.e. (a) first-order difference, and (b) half-wave rectification.\nSum across all frequency bins,  k .\nLuckily, librosa has librosa.onset.onset_strength which computes a novelty function using spectral flux.","metadata":{}},{"cell_type":"code","source":"spectral_novelty = librosa.onset.onset_strength(x, sr=sr)\nframes = numpy.arange(len(spectral_novelty))\nt = librosa.frames_to_time(frames, sr=sr)\nplt.figure(figsize=(15, 4))\nplt.plot(t, spectral_novelty, 'r-')\nplt.xlim(0, t.max())\nplt.xlabel('Time (sec)')\nplt.legend(('Spectral Novelty',))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Questions\nNovelty functions are dependent on frame_length and hop_length. Adjust these two parameters. How do they affect the novelty function?\n\nTry with other audio files. How do the novelty functions compare?","metadata":{}},{"cell_type":"markdown","source":"# Peak Picking\nPeak picking is the act of locating peaks in a signal. For example, in onset detection, we may want to find peaks in a novelty function. These peaks would correspond to the musical onsets.\n\nLet's load an example audio file.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nprint(x.shape, sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Listen to the audio file:\n\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the signal:\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Compute an onset envelope:\n\nhop_length = 256\nonset_envelope = librosa.onset.onset_strength(x, sr=sr, hop_length=hop_length)\nonset_envelope.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Generate a time variable:\n\nN = len(x)\nT = N/float(sr)\nt = numpy.linspace(0, T, len(onset_envelope))\n#Plot the onset envelope:\n\nplt.figure(figsize=(14, 5))\nplt.plot(t, onset_envelope)\nplt.xlabel('Time (sec)')\nplt.xlim(xmin=0)\nplt.ylim(0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"In this onset strength envelope, we clearly see many peaks. Some correspond to onsets, and others don't. How do we create peak picker that will detect true peaks while avoiding unwanted spurious peaks?\n\nlibrosa.util has a peak_pick method. We can control the parameters based upon our signal. Let's see how it works:","metadata":{}},{"cell_type":"markdown","source":"def peak_pick(x, pre_max, post_max, pre_avg, post_avg, delta, wait):\n    '''Uses a flexible heuristic to pick peaks in a signal.\n\n    A sample n is selected as a peak if the corresponding x[n]\n    fulfills the following three conditions:\n\n    1. `x[n] == max(x[n - pre_max:n + post_max])`\n    2. `x[n] >= mean(x[n - pre_avg:n + post_avg]) + delta`\n    3. `n - previous_n > wait`\n\n    where `previous_n` is the last sample picked as a peak (greedily).\nThis implementation is based on Boeck, Sebastian, Florian Krebs, and Markus Schedl. “Evaluating the Online Capabilities of Onset Detection Methods.” ISMIR. 2012.\n\nGet the frame indices of the peaks:","metadata":{}},{"cell_type":"code","source":"onset_frames = librosa.util.peak_pick(onset_envelope, 7, 7, 7, 7, 0.5, 5)\nonset_frames","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the onset envelope along with the detected peaks:\n\nplt.figure(figsize=(14, 5))\nplt.plot(t, onset_envelope)\nplt.grid(False)\nplt.vlines(t[onset_frames], 0, onset_envelope.max(), color='r', alpha=0.7)\nplt.xlabel('Time (sec)')\nplt.xlim(0, T)\nplt.ylim(0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#superimpose a click track upon the original:\n\nclicks = librosa.clicks(frames=onset_frames, sr=22050, hop_length=hop_length, length=N)\nipd.Audio(x+clicks, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Using the parameters above, we find that the peak picking algorithm seems to have high precision, e.g. few false positives. However, recall can be improved, i.e. it is missing several onsets that actually occur in the audio signal.","metadata":{}},{"cell_type":"markdown","source":"Questions\nAdjust the hop length from 512 to 256 or 1024. How does that affect the onset envelope, and consequently, the peak picking?\n\nAdjust the peak_pick parameters, pre_max, post_max, pre_avg, post_avg, delta, and wait. How do the detected peaks change?\n\nTry this notebook again on other audio files:","metadata":{}},{"cell_type":"markdown","source":"# Onset Detection\nAutomatic detection of musical events in an audio signal is one of the most fundamental tasks in music information retrieval. Here, we will show how to detect an onset, the very instant that marks the beginning of the transient part of a sound, or the earliest moment at which a transient can be reliably detected.\n\nFor more reading:\n\nTutorial on Onset Detection by Juan Bello\nBoeck, Sebastian, Florian Krebs, and Markus Schedl. “Evaluating the Online Capabilities of Onset Detection Methods.” ISMIR 2012\nLoad an audio file into the NumPy array x and sampling rate sr.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"librosa.onset.onset_detect works in the following way:\n\nCompute a spectral novelty function.\nFind peaks in the spectral novelty function.\n[optional] Backtrack from each peak to a preceding local minimum. Backtracking can be useful for finding segmentation points such that the onset occurs shortly after the beginning of the segment.\nCompute the frame indices for estimated onsets in a signal:","metadata":{}},{"cell_type":"code","source":"onset_frames = librosa.onset.onset_detect(x, sr=sr, wait=1, pre_avg=1, post_avg=1, pre_max=1, post_max=1)\nprint(onset_frames) # frame numbers of estimated onsets","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Convert onsets to units of seconds:\n\nonset_times = librosa.frames_to_time(onset_frames)\nprint(onset_times)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the onsets on top of a spectrogram of the audio:\n\nS = librosa.stft(x)\nlogS = librosa.amplitude_to_db(abs(S))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's also plot the onsets with the time-domain waveform.","metadata":{}},{"cell_type":"markdown","source":"# librosa.clicks\nWe can add a click at the location of each detected onset","metadata":{}},{"cell_type":"code","source":"clicks = librosa.clicks(frames=onset_frames, sr=sr, length=len(x))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Listen to the original audio plus the detected onsets. One way is to add the signals together, sample-wise:","metadata":{}},{"cell_type":"code","source":"ipd.Audio(x + clicks, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Another method is to play the original track in one stereo channel and the click track in the other stereo channel:","metadata":{}},{"cell_type":"code","source":"ipd.Audio(numpy.vstack([x, clicks]), rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#You can also change the click to a custom audio file instead:\n\ncowbell, _ = librosa.load('XC134874.mp3')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#More cowbell?\n\nclicks = librosa.clicks(frames=onset_frames, sr=sr, length=len(x), click=cowbell)\nipd.Audio(x + clicks, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Questions\nIn librosa.onset.onset_detect, use the backtrack=True parameter. What does that do, and how does it affect the detected onsets? (See librosa.onset.onset_backtrack.)\n\nIn librosa.onset.onset_detect, you can use the keyword parameters found in librosa.util.peak_pick, e.g. pre_max, post_max, pre_avg, post_avg, delta, and wait, to control the peak picking algorithm. Adjust these parameters. How does it affect the detected onsets?","metadata":{}},{"cell_type":"markdown","source":"# Onset-based Segmentation with Backtracking\nlibrosa.onset.onset_detect works by finding peaks in a spectral novelty function. However, these peaks may not actually coincide with the initial rise in energy or how we perceive the beginning of a musical note.\n\nThe optional keyword parameter backtrack=True will backtrack from each peak to a preceding local minimum. Backtracking can be useful for finding segmentation points such that the onset occurs shortly after the beginning of the segment. We will use backtrack=True to perform onset-based segmentation of a signal.\n\nLoad an audio file into the NumPy array x and sampling rate sr.","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nprint(x.shape, sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Listen:\n\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Compute the frame indices for estimated onsets in a signal:\nhop_length = 512\nonset_frames = librosa.onset.onset_detect(x, sr=sr, hop_length=hop_length)\nprint(onset_frames) \n# frame numbers of estimated onsets","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Convert onsets to units of seconds:\nonset_times = librosa.frames_to_time(onset_frames, sr=sr, hop_length=hop_length)\nprint(onset_times)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Convert onsets to units of samples:\nonset_samples = librosa.frames_to_samples(onset_frames, hop_length=hop_length)\nprint(onset_samples)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the onsets on top of a spectrogram of the audio:\nS = librosa.stft(x)\nlogS = librosa.logamplitude(S)\nlibrosa.display.specshow(logS, sr=sr, x_axis='time', y_axis='log')\nplt.vlines(onset_times, 0, 10000, color='k')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we see in the spectrogram, the detected onsets seem to occur a bit before the actual rise in energy.\n\nLet's listen to these segments. We will create a function to do the following:\n\nDivide the signal into segments beginning at each detected onset.\nPad each segment with 500 ms of silence.\nConcatenate the padded segments.","metadata":{}},{"cell_type":"code","source":"def concatenate_segments(x, onset_samples, pad_duration=0.500):\n    \"\"\"Concatenate segments into one signal.\"\"\"\n    silence = numpy.zeros(int(pad_duration*sr)) # silence\n    frame_sz = min(numpy.diff(onset_samples))   # every segment has uniform frame size\n    return numpy.concatenate([\n        numpy.concatenate([x[i:i+frame_sz], silence]) # pad segment with silence\n        for i in onset_samples\n    ])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Concatenate the segments:\nconcatenated_signal = concatenate_segments(x, onset_samples, 0.500)\n#Listen to the concatenated signal\nipd.Audio(concatenated_signal, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we hear, the little glitch between segments occurs because the segment boundaries occur during the attack, not before the attack.\n\nlibrosa.onset.onset_backtrack\nWe can avoid this glitch by backtracking from the detected onsets.\n\nWhen setting the parameter backtrack=True, librosa.onset.onset_detect will call librosa.onset.onset_backtrack. For each detected onset, librosa.onset.onset_backtrack searches backward for a local minimum.","metadata":{}},{"cell_type":"code","source":"onset_frames = librosa.onset.onset_detect(x, sr=sr, hop_length=hop_length, backtrack=True)\n#Convert onsets to units of seconds:\n\nonset_times = librosa.frames_to_time(onset_frames, sr=sr, hop_length=hop_length)\n#Convert onsets to units of samples:\n\nonset_samples = librosa.frames_to_samples(onset_frames, hop_length=hop_length)\n#Plot the onsets on top of a spectrogram of the audio:\n\nS = librosa.stft(x)\nlogS = librosa.logamplitude(S)\nlibrosa.display.specshow(logS, sr=sr, x_axis='time', y_axis='log')\nplt.vlines(onset_times, 0, 10000, color='k')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Notice how the vertical lines denoting each segment boundary appears before each rise in energy.\n\nConcatenate the segments:","metadata":{}},{"cell_type":"code","source":"concatenated_signal = concatenate_segments(x, onset_samples, 0.500)\n#Listen to the concatenated signal:\n\nipd.Audio(concatenated_signal, rate=sr)\n#While listening, notice now the segments are perfectly segmented.","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Tempo Estimation\nTempo (Wikipedia) refers to the speed of a musical piece. More precisely, tempo refers to the rate of the musical beat and is given by the reciprocal of the beat period. Tempo is often defined in units of beats per minute (BPM).\n\nIn classical music, common tempo markings include grave, largo, lento, adagio, andante, moderato, allegro, vivace, and presto. See Basic tempo markings for more.\n\n# Tempogram\nTempo can vary locally within a piece. Therefore, we introduce the tempogram (FMP, p. 317) as a feature matrix which indicates the prevalence of certain tempi at each moment in time.\n\n# Fourier Tempogram\nThe Fourier Tempogram (FMP, p. 319) is basically the magnitude spectrogram of the novelty function.\n\nLoad an audio file:","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The tempo of this excerpt is about 58/116 BPM.\n\nCompute the onset envelope, i.e. novelty function:","metadata":{}},{"cell_type":"code","source":"hop_length = 200 # samples per frame\nonset_env = librosa.onset.onset_strength(x, sr=sr, hop_length=hop_length, n_fft=2048)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the onset envelope:\n\nframes = range(len(onset_env))\nt = librosa.frames_to_time(frames, sr=sr, hop_length=hop_length)\nplt.plot(t, onset_env)\nplt.xlim(0, t.max())\nplt.ylim(0)\nplt.xlabel('Time (sec)')\nplt.title('Novelty Function')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Compute the short-time Fourier transform (STFT) of the novelty function. Since the novelty function is computed in frame increments, the hop length of this STFT should be pretty small:\n\nS = librosa.stft(onset_env, hop_length=1, n_fft=512)\nfourier_tempogram = numpy.absolute(S)\n#Plot the Fourier tempogram:\n\nlibrosa.display.specshow(fourier_tempogram, sr=sr, hop_length=hop_length, x_axis='time')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Autocorrelation Tempogram\nConsider a segment from the above novelty function:","metadata":{}},{"cell_type":"code","source":"n0 = 100\nn1 = 500\nplt.plot(t[n0:n1], onset_env[n0:n1])\nplt.xlim(t[n0], t[n1])\nplt.xlabel('Time (sec)')\nplt.title('Novelty Function')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the autocorrelation of this segment:\n\ntmp = numpy.log1p(onset_env[n0:n1])\nr = librosa.autocorrelate(tmp)\nplt.plot(t[:n1-n0], r)\nplt.xlim(t[0], t[n1-n0])\nplt.xlabel('Lag (sec)')\nplt.ylim(0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Wherever the autocorrelation is high is a good candidate of the beat period.\n\nplt.plot(60/t[:n1-n0], r)\nplt.xlim(20, 200)\nplt.xlabel('Tempo (BPM)')\nplt.ylim(0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#We will apply this principle of autocorrelation to estimate the tempo at every segment in the novelty function.\n\n#librosa.feature.tempogram implements an autocorrelation tempogram, a short-time autocorrelation of the (spectral) novelty function.\n\n#For more information:\n\n#Grosche, Peter, Meinard Müller, and Frank Kurth. “Cyclic tempogram - A mid-level tempo representation for music signals.” ICASSP, 2010.\n#Compute a tempogram:\n\ntempogram = librosa.feature.tempogram(onset_envelope=onset_env, sr=sr, hop_length=hop_length, win_length=400)\nlibrosa.display.specshow(tempogram, sr=sr, hop_length=hop_length, x_axis='time', y_axis='tempo')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Estimating Global Tempo¶\nWe will use librosa.beat.tempo to estimate the global tempo in an audio file.\n\nEstimate the tempo:","metadata":{}},{"cell_type":"code","source":"tempo = librosa.beat.tempo(x, sr=sr)\nprint(tempo)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Visualize the tempo estimate on top of the input signal:\n\nT = len(x)/float(sr)\nseconds_per_beat = 60.0/tempo[0]\nbeat_times = numpy.arange(0, T, seconds_per_beat)\nlibrosa.display.waveplot(x)\nplt.vlines(beat_times, -1, 1, color='r')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Listen to the input signal with a click track using the tempo estimate:\n\nclicks = librosa.clicks(beat_times, sr, length=len(x))\nipd.Audio(x + clicks, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# BEAT Tracking\nlibrosa.beat.beat_track\nLoad an audio file:","metadata":{}},{"cell_type":"code","source":"x, sr = librosa.load('XC134874.mp3')\nipd.Audio(x, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Use librosa.beat.beat_track to estimate the beat locations and the global tempo:","metadata":{}},{"cell_type":"code","source":"tempo, beat_times = librosa.beat.beat_track(x, sr=sr, start_bpm=60, units='time')\nprint(tempo)\nprint(beat_times)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot the beat locations over the waveform:\n\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, alpha=0.6)\nplt.vlines(beat_times, -1, 1, color='r')\nplt.ylim(-1, 1)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Plot a histogram of the intervals between adjacent beats:\n\nbeat_times_diff = numpy.diff(beat_times)\nplt.figure(figsize=(14, 5))\nplt.hist(beat_times_diff, bins=50, range=(0,4))\nplt.xlabel('Beat Length (seconds)')\nplt.ylabel('Count')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Visually, it's difficult to tell how correct the estimated beats are. Let's listen to a click track:\n\nclicks = librosa.clicks(beat_times, sr=sr, length=len(x))\nipd.Audio(x + clicks, rate=sr)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Use the IPython interactive widgets to observe how the output changes as we vary the parameters of the beat tracker.","metadata":{}},{"cell_type":"code","source":"def f(start_bpm, tightness_exp):\n    return librosa.beat.beat_track(x, sr=sr, start_bpm=start_bpm, tightness=10**tightness_exp, units='time')\ninteract(f, start_bpm=60, tightness_exp=2)","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}