Quantify amount of shaking caused by motorized liquid injection during pTIRF imaging

quantify_shaking.py

#!/usr/bin/env python
# coding: utf-8

import skimage.io as skio
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from os.path import basename, splitext, exists
from time import sleep
from os import listdir


def quantify_shaking(tiff_file, crop=10):
    # read tiff file
    img = skio.imread(tiff_file)

    # crop the central part
    N, h, w = img.shape
    _from, _to = 1 - crop / 100.0, 1 + crop / 100.0
    img = img[
        :, int(h / 2 * _from) : int(h / 2 * _to), int(w / 2 * _from) : int(w / 2 * _to)
    ]
    Npx = np.prod(img[0].shape)

    # compute image intensity
    img_I = img.mean(axis=(1, 2))
    img_I = img_I - img_I.min() + img_I.ptp() / 20

    # compute normalized intraframe variance
    img_intraframe_v = img.std(axis=(1, 2)) ** 2 / Npx / img_I

    # compute interframe variance
    diff = np.diff(img.astype("float"), axis=0)
    img_interframe_v = np.std(diff, axis=(1, 2)) ** 2 / Npx

    # compute ratio of interframe to intraframe variance, then
    # remove the baseline by taking a derivative and calculating a square
    return img_I, np.gradient(img_interframe_v / img_intraframe_v[1:]) ** 2


def plot_shaking(tiff_file, **kwargs):
    img_I, shaking = quantify_shaking(tiff_file, **kwargs)
    fname = splitext(basename(tiff_file))[0]
    fig, ax1 = plt.subplots()
    ax2 = ax1.twinx()
    ax2.plot(img_I, linewidth=0.75, alpha=0.6, color="darkred")
    ax2.set_ylabel("Intensity")

    ax1.plot(shaking, linewidth=0.75)
    plt.title(
        f"{fname}\nShaking intensity: {sum(shaking):.1f} a.u.  (in {len(np.where(shaking > 10)[0])} frames)"
    )
    plt.xlabel("Frame number")
    ax1.set_ylabel("Shaking, a.u.")
    fig.patch.set_alpha(1.0)
    plt.savefig(fname + "_shaking.png", dpi=600, transparent=False)
    plt.close()

    pd.DataFrame(dict(
        intensity=img_I[1:],
        shaking=shaking)
    ).to_csv(fname + "_shaking.csv")


def main():
    while True:
        sleep(0.1)
        try:
            for f in listdir("."):
                name, ext = splitext(f)
                if ".tif" in ext.lower():
                    if not exists(f"{name}_shaking.csv"):
                        plot_shaking(f)
                        print(f)
        except PermissionError:
            pass
        except KeyboardInterrupt:
            break


if __name__ == "__main__":
    main()