main.py

#!/usr/bin/env python

"""
MPAI CAE-ARP Tape Audio Restoration.

Implements MPAI CAE-ARP Tape Audio Restoration Technical Specification.
It identifies and restore portions of the Preservation Audio File, providing:
- Restored Audio Files;
- Editing List
"""

import matplotlib.pyplot as plt
import numpy as np
import os
import shutil
import yaml
from control import c2d, TransferFunction
from scipy.io import wavfile
from scipy.signal import freqs, freqz, tf2zpk, zpk2tf, lfilter

__author__ = "Nadir Dalla Pozza"
__copyright__ = "Copyright 2022, Audio Innova S.r.l."
__credits__ = ["Niccolò Pretto", "Nadir Dalla Pozza", "Sergio Canazza"]
__license__ = "GPL v3.0"
__version__ = "1.0.1"
__maintainer__ = "Nadir Dalla Pozza"
__email__ = "nadir.dallapozza@unipd.it"
__status__ = "Production"


# Class for customizing Console Colors
class CC:
    PURPLE = '\033[95m'
    CYAN = '\033[96m'
    DARK_CYAN = '\033[36m'
    BLUE = '\033[94m'
    GREEN = '\033[92m'
    YELLOW = '\033[93m'
    RED = '\033[91m'
    BOLD = '\033[1m'
    UNDERLINE = '\033[4m'
    END = '\033[0m'


def save_file(file):
    if not os.path.exists(temp_path):
        # Create directory
        os.mkdir(temp_path)
        print("temp directory '% s' created" % temp_path)
    raf_path = os.path.join(temp_path, 'RestoredAudioFiles')
    make_raf = False
    if not os.path.exists(raf_path):
        # Create directory
        os.mkdir(raf_path)
        make_raf = True
        print("Restored Audio Files directory '% s' created" % raf_path)
    else:
        print((CC.PURPLE + "Restored Audio Files directory '% s' already exists!" + CC.END) % raf_path)
        overwrite = input('Do you want to overwrite it? [y/n]: ')
        if overwrite.casefold() == 'y':
            # Overwrite directory
            shutil.rmtree(raf_path)
            os.mkdir(raf_path)
            make_raf = True
            print('Restored Audio Files directory overwritten')
        elif overwrite.casefold() != 'n':
            print(CC.RED + 'Unknown command, exiting' + CC.END)
            quit(os.EX_USAGE)
    if make_raf:
        print("Saving Restored Audio File to: '%s' ..." % raf_path)
        wavfile.write(os.path.join(raf_path, '1.wav'), FS, file)


if __name__ == '__main__':
    # Read configuration file
    config = object
    try:
        config = yaml.safe_load(open('config.yaml', 'r'))
        if 'WORKING_PATH' not in config:
            print(CC.RED + 'WORKING_PATH key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'PRESERVATION_FILE_NAME' not in config:
            print(CC.RED + 'PRESERVATION_FILE_NAME key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'STANDARD_W' not in config:
            print(CC.RED + 'STANDARD_W key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'SPEED_W' not in config:
            print(CC.RED + 'SPEED_W key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'STANDARD_R' not in config:
            print(CC.RED + 'STANDARD_R key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'SPEED_R' not in config:
            print(CC.RED + 'SPEED_R key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'FS' not in config:
            print(CC.RED + 'FS key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'BPS' not in config:
            print(CC.RED + 'BPS key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
        if 'PLOTS' not in config:
            print(CC.RED + 'PLOTS key not found in config.yaml!' + CC.END)
            quit(os.EX_CONFIG)
    except FileNotFoundError:
        print(CC.RED + 'config.yaml file not found!' + CC.END)
        quit(os.EX_NOINPUT)

    # config.yaml variables
    WORKING_PATH = config['WORKING_PATH']
    PRESERVATION_FILE_NAME = config['PRESERVATION_FILE_NAME']
    STANDARD_W = config['STANDARD_W']
    SPEED_W = config['SPEED_W']
    STANDARD_R = config['STANDARD_R']
    SPEED_R = config['SPEED_R']
    FS = config['FS']
    BPS = config['BPS']
    PLOTS = config['PLOTS']

    # Output path
    temp_path = os.path.join(WORKING_PATH, 'temp', PRESERVATION_FILE_NAME)

    # Configuration parameters check

    # Recording tape speed check
    if SPEED_W != 3.75 and SPEED_W != 7.5 and SPEED_W != 15 and SPEED_W != 30:
        print(
            CC.RED + 'Incorrect SPEED_W: \'' + str(SPEED_W) + '\'. Accepted value are: 3.75, 7.5, 15, 30.' + CC.END
        )
        quit(os.EX_CONFIG)

    # Reading tape speed check.
    if SPEED_R != 3.75 and SPEED_R != 7.5 and SPEED_R != 15 and SPEED_R != 30:
        print(
            CC.RED + 'Incorrect SPEED_R: \'' + str(SPEED_R) + '\'. Accepted value are: 3.75, 7.5, 15, 30.' + CC.END
        )
        quit(os.EX_CONFIG)

    # Equalization standard check.
    if not (STANDARD_R == 'CCIR' or STANDARD_R == 'NAB'):
        print(
            CC.RED + 'Incorrect STANDARD_R: \'' + STANDARD_R + '\'. Accepted values are: CCIR, NAB.' + CC.END
        )
        quit(os.EX_CONFIG)
    if not (STANDARD_W == 'CCIR' or STANDARD_W == 'NAB'):
        print(
            CC.RED + 'Incorrect STANDARD_W: \'' + STANDARD_W + '\'. Accepted values are: CCIR, NAB.' + CC.END
        )
        quit(os.EX_CONFIG)

    # CCIR speed check.
    if STANDARD_W == 'CCIR' and SPEED_W == 3.75:
        print(
            CC.YELLOW + 'CCIR is undefined at 3.75 ips. Recording equalization standard is set to NAB.' + CC.END
        )
        STANDARD_W = 'NAB'
    if STANDARD_R == 'CCIR' and SPEED_R == 3.75:
        print(
            CC.YELLOW + 'CCIR is undefined at 3.75 ips. Reading equalization standard is set to NAB.' + CC.END
        )
        STANDARD_R = 'NAB'
    # NAB speed check.
    if STANDARD_W == 'NAB' and SPEED_W == 30:
        print(
            CC.YELLOW + 'NAB is undefined at 30 ips. Recording equalization standard is set to CCIR.' + CC.END
        )
        STANDARD_W = 'CCIR'
    if STANDARD_R == 'NAB' and SPEED_R == 30:
        print(
            CC.YELLOW + 'NAB is undefined at 30 ips. Reading equalization standard is set to CCIR.' + CC.END
        )
        STANDARD_R = 'CCIR'

    # Sampling frequency check.
    if len(PRESERVATION_FILE_NAME) == 0:
        if FS != 44100 and FS != 48000 and FS != 96000:
            print(
                CC.RED + 'Incorrect FS: \'' + str(FS) + '\'. Accepted values are: 44100, 48000, 96000.' + CC.END
            )
            quit(os.EX_CONFIG)

    # Bits per sample check.
    if BPS != 8 and BPS != 16 and BPS != 24 and BPS != 32:
        print(
            CC.RED + 'Incorrect BPS: \'' + str(BPS) + '\'. Accepted values are: 8, 16, 24, 32.' + CC.END
        )
        quit(os.EX_CONFIG)

    # Display input parameters
    print('Input parameters:')
    print('    WORKING_PATH:           ' + WORKING_PATH)
    print('    PRESERVATION_FILE_NAME: ' + PRESERVATION_FILE_NAME)
    print('    STANDARD_W:             ' + STANDARD_W)
    print('    SPEED_W:                ' + str(SPEED_W) + ' ips')
    print('    STANDARD_R:             ' + STANDARD_R)
    print('    SPEED_R:                ' + str(SPEED_R) + ' ips')
    if len(PRESERVATION_FILE_NAME) == 0:
        print('    FS:                     ' + str(FS) + ' Hz')
    print('    BPS:                    ' + str(BPS) + '\n')

    # Preservation Audio File check
    paf = []
    if len(PRESERVATION_FILE_NAME) > 0:
        audio_file = PRESERVATION_FILE_NAME + '.wav'
        paf_path = os.path.join(WORKING_PATH, 'PreservationAudioFile', audio_file)
        try:
            print("Opening '%s'..." % paf_path)
            FS, paf = wavfile.read(paf_path)
            print('Preservation Audio File opened!')
            print('Overwritten parameters:')
            print('    FS:                     ' + str(FS) + ' Hz\n')
        except OSError:
            print(CC.RED + "Preservation Audio File not found!" + CC.END)
            quit(os.EX_NOINPUT)

    # Equalization standard time constants

    # CCIR time constants.
    t2_30 = 17.5 * 10**(-6)  # time constant CCIR_30
    t2_15 = 35 * 10**(-6)  # time constant CCIR_15
    t2_7 = 70 * 10**(-6)  # time constant CCIR_7.5

    # NAB time constants.
    t3 = 3180 * 10**(-6)
    t4_15 = 50 * 10**(-6)  # time constant NAB_15
    t4_7 = 50 * 10**(-6)  # time constant NAB_7.5
    t4_3 = 90 * 10**(-6)  # time constant NAB_3.75

    # Decision stage

    a = []  # Filter numerator
    b = []  # Filter denominator
    case = 0  # Reference case number

    # This section will establish which time constants must be modified to obtain the desired equalisation standard.
    if STANDARD_W == 'CCIR':
        if SPEED_W == 30:
            if STANDARD_R == 'NAB':
                # Case 1
                if SPEED_R == 15:
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Correction filter: NABw15_mod + CCIRr30
                    # - NAB constants divided by 2
                    t3 = t3 / 2
                    t4 = t4_15 / 2
                    # - CCIR_30 constant not altered
                    t2 = t2_30
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 1
                # Case 2
                elif SPEED_R == 7.5:
                    FS = 4 * FS  # Quadrupling the sampling frequency
                    # Correction filter: NABw7.5_mod + CCIRr30
                    # - NAB constants divided by 4
                    t3 = t3 / 4
                    t4 = t4_7 / 4
                    # - CCIR_30 constant not altered
                    t2 = t2_30
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 2
                # Case 3
                else:  # SPEED_R == 3.75
                    FS = 8 * FS  # Multiplying by 8 the sampling frequency
                    # Correction filter: NABw3.75_mod + CCIRr30
                    # - NAB constants divided by 8
                    t3 = t3 / 8
                    t4 = t4_3 / 8
                    # - CCIR_30 constant not altered
                    t2 = t2_30
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 3
            else:  # STANDARD_R == 'CCIR'
                # Case 31
                if SPEED_R == 30:
                    print('Reference case: 31')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
                # Case 15
                elif SPEED_R == 15:
                    FS = 2 * FS  # Doubling sampling frequency
                    # Plot information
                    case = 15
                # Case 16
                else:  # SPEED_R == 7.5
                    FS = 4 * FS  # Quadrupling the sampling frequency
                    # Plot information
                    case = 16
        elif SPEED_W == 15:
            if STANDARD_R == 'NAB':
                # Case 28
                if SPEED_R == 15:
                    # No speed change
                    # Correction filter: NABw15 + CCIRr15
                    # - NAB_15 constants not altered
                    t4 = t4_15
                    # - CCIR_30 constant not altered
                    t2 = t2_15
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 28
                # Case 6
                elif SPEED_R == 7.5:
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Correction filter: NABw7.5_mod + CCIRr15
                    # - NAB constants divided by 2
                    t3 = t3 / 2
                    t4 = t4_7 / 2
                    # - CCIR_15 constant not altered
                    t2 = t2_15
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 6
                # Case 7
                else:  # SPEED_R == 3.75
                    FS = 4 * FS  # Quadrupling the sampling frequency
                    # Correction filter: NABw3.75_mod + CCIRr15
                    # - NAB constants divided by 4
                    t3 = t3 / 4
                    t4 = t4_3 / 4
                    # - CCIR_15 constant not altered
                    t2 = t2_15
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 7
            else:  # STANDARD_R == 'CCIR'
                # Case 19
                if SPEED_R == 30:
                    FS = FS / 2  # Halving the sampling frequency
                    # Plot information
                    case = 19
                # Case 33
                elif SPEED_R == 15:
                    print('Reference case: 33')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
                # Case 20
                else:  # SPEED_R == 7.5
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Plot information
                    case = 20
        else:  # SPEED_W == 7.5
            if STANDARD_R == 'NAB':
                # Case 10
                if SPEED_R == 15:
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: NABw15_mod + CCIRr7.5
                    # - NAB constants multiplied by 2
                    t3 = t3 * 2
                    t4 = t4_15 * 2
                    # - CCIR_7.5 constant not altered
                    t2 = t2_7
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 10
                # Case 30
                elif SPEED_R == 7.5:
                    # No speed change
                    # Correction filter: NABw7.5 + CCIRr7.5
                    # - NAB_7.5 constant not altered
                    t4 = t4_7
                    # - CCIR_7.5 constant not altered
                    t2 = t2_7
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 30
                # Case 11
                else:  # SPEED_R == 3.75
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Correction filter: NABw3.75_mod + CCIRr7.5
                    # - NAB constants divided by 2
                    t3 = t3 / 2
                    t4 = t4_3 / 2
                    # - CCIR_7.5 constant not altered
                    t2 = t2_7
                    # Filter coefficients
                    a = [t2 * t3, t2 + t3, 1]
                    b = [t3 * t4, t3, 0]
                    # Plot information
                    case = 11
            else:  # STANDARD_R == 'CCIR'
                # Case 23
                if SPEED_R == 30:
                    FS = FS / 4  # Quartering the sampling frequency
                    # Plot information
                    case = 23
                # Case 24
                elif SPEED_R == 15:
                    FS = FS / 2  # Halving the sampling frequency
                    # Plot information
                    case = 24
                # Case 35
                else:  # SPEED_R == 7.5
                    print('Reference case: 35')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
    else:  # STANDARD_W == 'NAB'
        if SPEED_W == 15:
            if STANDARD_R == 'NAB':
                # Case 32
                if SPEED_R == 15:
                    print('Reference case: 32')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
                # Case 17
                elif SPEED_R == 7.5:
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Correction filter: NABw7.5_mod + NABr15
                    # - NABw constants divided by 2
                    t3_mod = t3 / 2
                    t4_mod = t4_7 / 2
                    # - NABr constant not altered
                    t4 = t4_15
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 17
                # Case 18
                else:  # SPEED_R == 3.75
                    FS = 4 * FS  # Quadrupling the sampling frequency
                    # Correction filter: NABw3.75_mod + NABr15
                    # - NAB constants divided by 4
                    t3_mod = t3 / 4
                    t4_mod = t4_3 / 4
                    # - NABr constant not altered
                    t4 = t4_15
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 18
            else:  # STANDARD_R == 'CCIR'
                # Case 4
                if SPEED_R == 30:
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: CCIRw30_mod + NABr15
                    # - CCIR_30 constant multiplied by 2
                    t2 = t2_30 * 2
                    # - NAB_15 constant not altered
                    t4 = t4_15
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 4
                # Case 27
                elif SPEED_R == 15:
                    # No speed change
                    # Correction filter: CCIRw15 + NABr15
                    # - CCIR_15 constant not altered
                    t2 = t2_15
                    # - NAB_15 constant not altered
                    t4 = t4_15
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 27
                # Case 5
                else:  # SPEED_R == 7.5
                    FS = FS * 2  # Doubling the sampling frequency
                    # Correction filter: CCIRw7.5_mod + NABr15
                    # - CCIR_7.5 constant divided by 2
                    t2 = t2_7 / 2
                    # - NAB_15 constant not altered
                    t4 = t4_15
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 5
        elif SPEED_W == 7.5:
            if STANDARD_R == 'NAB':
                # Case 21
                if SPEED_R == 15:
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: NABw15_mod + NABr7.5
                    # - NABw constants multiplied by 2
                    t3_mod = t3 * 2
                    t4_mod = t4_15 * 2
                    # - NABr constant not altered
                    t4 = t4_7
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 21
                # Case 34
                elif SPEED_R == 7.5:
                    print('Reference case: 34')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
                # Case 22
                else:  # SPEED_R == 3.75
                    FS = 2 * FS  # Doubling the sampling frequency
                    # Correction filter: NABw3.75_mod + NABr7.5
                    # - NABw constants divided by 2
                    t3_mod = t3 / 2
                    t4_mod = t4_3 / 2
                    # - NABr constant not altered
                    t4 = t4_7
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 22
            else:  # STANDARD_R == 'CCIR'
                # Case 8
                if SPEED_R == 30:
                    FS = FS / 4  # Quartering the sampling frequency
                    # Correction filter: CCIRw30_mod + NABr7.5
                    # - CCIR_30 constant multiplied by 4
                    t2 = t2_30 * 4
                    # - NAB_7.5 constant not altered
                    t4 = t4_7
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 8
                # Case 9
                elif SPEED_R == 15:
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: CCIRw15_mod + NABr7.5
                    # - CCIR_15 constant multiplied by 2
                    t2 = t2_15 * 2
                    # - NAB_7.5 constant not altered
                    t4 = t4_7
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 9
                # Case 29
                else:  # SPEED_R == 7.5
                    # No speed change
                    # Correction filter: CCIRw7.5 + NABr7.5
                    # - CCIR_7.5 constant not altered
                    t2 = t2_7
                    # - NAB_7.5 constant not altered
                    t4 = t4_7
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 29
        else:  # SPEED_W == 3.75
            if STANDARD_R == 'NAB':
                # Case 25
                if SPEED_R == 15:
                    FS = FS / 4  # Quartering the sampling frequency
                    # Correction filter: NABw15_mod + NABr3.75
                    # - NAB constants multiplied by 4
                    t3_mod = t3 * 4
                    t4_mod = t4_15 * 4
                    # - NABr constant not altered
                    t4 = t4_3
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 25
                # Case 26
                elif SPEED_R == 7.5:
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: NABw7.5_mod + NABr3.75
                    # - NAB constants multiplied by 2
                    t3_mod = t3 * 2
                    t4_mod = t4_7 * 2
                    # - NABr constant not altered
                    t4 = t4_3
                    # Filter coefficients
                    a = [t3 * t3_mod * t4, t3 * (t3_mod + t4), t3]
                    b = [t3 * t3_mod * t4_mod, t3_mod * (t3 + t4_mod), t3_mod]
                    # Plot information
                    case = 26
                # Case 36
                else:  # SPEED_R == 3.75
                    print('Reference case: 36')
                    print(CC.GREEN + 'Nothing to do!' + CC.END)
                    quit(os.EX_OK)
            else:  # STANDARD_R == 'CCIR'
                # Case 12
                if SPEED_R == 30:
                    FS = FS / 8  # Dividing by 8 the sampling frequency
                    # Correction filter: CCIRw30_mod + NABr3.75
                    # - CCIR_30 constant multiplied by 8
                    t2 = t2_30 * 8
                    # - NAB_3.75 constant not altered
                    t4 = t4_3
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 12
                # Case 13
                elif SPEED_R == 15:
                    FS = FS / 4  # Quartering the sampling frequency
                    # Correction filter: CCIRw15_mod + NABr3.75
                    # - CCIR_15 constant multiplied by 4
                    t2 = t2_15 * 4
                    # - NAB_3.75 constant not altered
                    t4 = t4_3
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 13
                # Case 14
                else:  # SPEED_R == 7.5
                    FS = FS / 2  # Halving the sampling frequency
                    # Correction filter: CCIRw7.5_mod + NABr3.75
                    # - CCIR_7.5 constant multiplied by 2
                    t2 = t2_7 * 2
                    # - NAB_3.75 constant not altered
                    t4 = t4_3
                    # Filter coefficients
                    a = [t3 * t4, t3, 0]
                    b = [t2 * t3, t2 + t3, 1]
                    # Plot information
                    case = 14

    # Casting FS to int because wavfile.write() is stupid
    FS = round(FS)
    print('Reference case: ' + str(case))
    print('Operational sampling frequency: ' + str(FS) + ' Hz.')

    # Correction filter

    # Not all cases present a correction filter!
    if len(a) != 0:
        # Analog transfer function
        H_a = TransferFunction(a, b)
        # Analog frequency vector
        w_a = np.logspace(np.log10(1), np.log10(FS * np.pi), 5000)

        if PLOTS:
            # Analog filter frequency response
            w_t, h_t = freqs(a, b, worN=w_a)
            # Plot analog graph
            # - Magnitude
            plt.subplot(2, 1, 1)
            plt.semilogx(w_t / (2 * np.pi), 20 * np.log10(abs(h_t)))
            plt.xlim([1, 24000])
            plt.xlabel('Frequency')
            plt.ylim([-40, 40])
            plt.ylabel('Amplitude response [dB]')
            plt.grid(True)
            # - Phase
            plt.subplot(2, 1, 2)
            plt.semilogx(w_t / (2 * np.pi), np.angle(h_t) * 180 / np.pi)
            plt.xlim([1, 24000])
            plt.xlabel('Frequency')
            plt.ylabel('Phase [deg]')
            plt.grid(True)

        # Digital transfer function through bilinear digitisation
        H_d = c2d(H_a, 1/FS, 'bilinear')
        num_d = H_d.num[0][0]  # Inspect Hd.num to see why [0][0] is needed...
        den_d = H_d.den[0][0]  # Same story here
        # Digital frequency vector
        w_d = np.logspace(np.log10(1), np.log10(FS / 2), 5000)

        if PLOTS:
            # Digital filter frequency response
            w_n, h_n = freqz(num_d, den_d, worN=w_d, fs=FS)
            # Plot digital graph
            # - Magnitude
            plt.subplot(2, 1, 1)
            plt.semilogx(w_n, 20 * np.log10(abs(h_n)), '--')
            plt.legend(['Analog', 'Bilinear'])
            # - Phase
            plt.subplot(2, 1, 2)
            plt.semilogx(w_n, np.angle(h_n) * 180 / np.pi, '--')
            plt.legend(['Analog', 'Bilinear'])

        # Pole check

        # New pole frequency
        pole_frequency = 2
        # Move to zero-pole representation
        z, p, k = tf2zpk(a, b)
        # Check if the function presents a pole at 0 Hz
        for i in range(len(p)):
            if p[i] == 0:
                # Replace pole
                p[i] = -pole_frequency * 2 * np.pi
                print('\n' + CC.PURPLE + 'Pole at 0 Hz replaced!' + CC.END)
                # Back to transfer function representation
                ap, bp = zpk2tf(z, p, k)

                # Analog transfer function
                Hp_a = TransferFunction(ap, bp)

                if PLOTS:
                    # Analog filter frequency response
                    wp_t, hp_t = freqs(ap, bp, worN=w_a)
                    # Plot analog graph
                    # - Magnitude
                    plt.subplot(2, 1, 1)
                    plt.semilogx(wp_t / (2 * np.pi), 20 * np.log10(abs(hp_t)))
                    # - Phase
                    plt.subplot(2, 1, 2)
                    plt.semilogx(wp_t / (2 * np.pi), np.angle(hp_t) * 180 / np.pi)

                # Digital transfer function through bilinear digitisation
                Hp_d = c2d(Hp_a, 1 / FS, 'bilinear')
                num_d = Hp_d.num[0][0]
                den_d = Hp_d.den[0][0]

                if PLOTS:
                    # Digital filter frequency response
                    wp_n, hp_n = freqz(num_d, den_d, worN=w_d, fs=FS)
                    # Plot digital graph
                    # - Magnitude
                    plt.subplot(2, 1, 1)
                    plt.semilogx(wp_n, 20 * np.log10(abs(hp_n)), '--')
                    plt.legend(['Analog', 'Bilinear', 'Pole - Analog', 'Pole - Digital'])
                    # - Phase
                    plt.subplot(2, 1, 2)
                    plt.semilogx(wp_n, np.angle(hp_n) * 180 / np.pi, '--')
                    plt.legend(['Analog', 'Bilinear', 'Pole - Analog', 'Pole - Digital'])

        if PLOTS:
            plt.show()

        if len(PRESERVATION_FILE_NAME) > 0:
            print('\nFiltering Preservation Audio File...')
            # Filter Preservation Audio File
            raf = lfilter(num_d, den_d, paf, axis=0)
            # Again, wavfile.write() is stupid, and you must cast everything to not destroy your ears...
            raf = np.rint(raf).astype(paf.dtype)

            # Save Restored Audio File
            save_file(raf)

    else:
        # Save Restored Audio File
        save_file(paf)

    # End
    print(CC.GREEN + CC.BOLD + "Success!" + CC.END + '\n')