raab/competition_experiments
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

all files for init commit added

Browse Source
This commit is contained in:
tillraab 2022-06-02 13:51:21 +02:00
commit a70ba50037
5 changed files with 563 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

81
LED_detect.py Normal file
View File

@ -0,0 +1,81 @@
import os
import cv2
import argparse
import numpy as np
import matplotlib.pyplot as plt
def check_LED(cap, frame_count, x0, x1, y0, y1):
fig, ax = plt.subplots()
ax.plot([x0, x0], [y0, y1], 'r')
ax.plot([x1, x1], [y0, y1], 'r')
ax.plot([x0, x1], [y0, y0], 'r')
ax.plot([x0, x1], [y1, y1], 'r')
plt.ion()
cap.set(cv2.CAP_PROP_POS_FRAMES, int(frame_count / 2))
f = None
try:
for i in np.arange(int(frame_count / 2), frame_count):
ret, frame = cap.read()
if f == None:
f = ax.imshow(frame)
else:
f.set_data(frame)
sum_frame = np.sum(frame, axis=2)
LED_v = np.mean(sum_frame[y0:y1, x0:x1])
print('%.0f: %.1f \n' % (i, LED_v))
plt.pause(0.001)
except KeyboardInterrupt:
plt.close()
quit()
quit()
def main(file_path, check, x, y):
folder, filename = os.path.split(file_path)
cap = cv2.VideoCapture(file_path)
frame_count = int(cap.get(cv2.CAP_PROP_FRAME_COUNT))
x0, x1 = x
y0, y1 = y
if check:
check_LED(cap, frame_count, x0, x1, y0, y1)
###############
# cap.set(cv2.CAP_PROP_POS_FRAMES, int(frame_count / 2))
# frame_count = 1000
########################
light_th = 100
LED_val = np.zeros(frame_count)
print('Frame_count: %.0f' % frame_count)
for i in range(frame_count):
if i % 1000 == 0:
print('progress: %.1f' % ((i/frame_count)*100) + '%')
ret, frame = cap.read()
sum_frame = np.sum(frame, axis=2)
LED_val[i] = np.mean(sum_frame[y0:y1, x0:x1])
np.save(os.path.join(folder, 'LED_val.npy'), LED_val)
LED_frames = np.arange(len(LED_val)-1)[(LED_val[:-1] < light_th) & (LED_val[1:] > light_th)]
np.save(os.path.join(folder, 'LED_frames.npy'), LED_frames)
fig, ax = plt.subplots()
ax.plot(np.arange(len(LED_val)), LED_val, color='k')
ax.plot(LED_frames, np.ones(len(LED_frames))*light_th, 'o', color='firebrick')
plt.show()
if __name__ == '__main__':
parser = argparse.ArgumentParser(description='Detect frames of blinking LED in video recordings.')
parser.add_argument('file', type=str, help='video file to be analyzed')
parser.add_argument("-c", action="store_true", help="check if LED pos is correct")
parser.add_argument('-x', type=int, nargs=2, default=[1272, 1282], help='x-borders of LED detect area (in pixels)')
parser.add_argument('-y', type=int, nargs=2, default=[1500, 1516], help='y-borders of LED area (in pixels)')
args = parser.parse_args()
import glob
main(args.file, args.c, args.x, args.y)

53
README.md Normal file
View File

@ -0,0 +1,53 @@
# How to competition experiment
__Workflow__ (Python-scripts/applications):
1) wavetracker.trackingGUI
2) wavetracker.EODsorter
3) LED_detect.py
4) eval_LED.py
5) trial_analysis.py
6) event_videos.py (optional)
## Raw data analysis using the wavetracker-modul
### trackingGUI.py
__Frequency extraction and tracking__
- open Raw-file (traces-grid1.raw)
- 'Spectrogram'-settings:
- overlap fraction: 0.8
- frequency resolution: 1
- check 'Auto-save'; press 'Run'
__Fine spectrogram__
- repeat steps above but press 'Calc. fine spec' instead of Run
- fine spec data saved in /home/"user"/analysis/"filename"
### EODsorter.py
- load dataset/folder
- correct tracked EOD traces
- fill EOD traces
- fine spec data needs to be manually added to the dataset-folder
## Competition trial analysis
### trail_analysis.py
- Detection of winners, their EODf traces, rises, etc. Results stored in "data-path"/analysis.
- (optional) Meta.csv file in base-path of analyzed data. Creates entries for each
analyzed recording (index = file names) and stores Meta-data. Manual competation suggested.
## Video analysis
### LED_detect.py
- Detect blinking LED (emitted by electric recording setup). Used for synchronization.
- "-c" argument to identify correct detection area for LED
- '-x' (tuple) borders of LED detection window on X-axis (in pixels)
- '-y' (tuple) borders of LED detection window on Y-axis (in pixels)
### eval_LED.py
- creates time vector to synchronize electric and video recording
- for each frame contains a time-point (in s) that corresponds to the electric recordings.
## Rise videos (optional)
- generates for each detected rise a short video showing the fish's behavior around the rise event.
- sorted in 'base-path'/rise_video.

59
eval_LED.py Normal file
View File

@ -0,0 +1,59 @@
import pandas as pd
import numpy as np
import sys
import os
import matplotlib.pyplot as plt
from IPython import embed
import cv2
import glob
def main(folder):
sr = 20000
video_path = glob.glob(os.path.join(folder, '2022*.mp4'))[0]
cap = cv2.VideoCapture(video_path)
frame_count = int(cap.get(cv2.CAP_PROP_FRAME_COUNT))
times = np.load(os.path.join(folder, 'times.npy'))
LED_idx = pd.read_csv(os.path.join(folder, 'led_idxs.csv'), sep=',')
led_idx = np.array(LED_idx).T[0]
led_frame = np.load(os.path.join(folder, 'LED_frames.npy'))
led_vals = np.load(os.path.join(folder, 'LED_val.npy'))
led_idx_span = led_idx[-1] - led_idx[0]
led_frame_span = led_frame[-1] - led_frame[0]
led_frame_to_idx = ((led_frame-led_frame[0]) / led_frame_span) * led_idx_span + led_idx[0]
frame_idxs = np.arange(frame_count)
frame_times = (((frame_idxs - led_frame[0]) / led_frame_span) * led_idx_span + led_idx[0]) / sr
if not os.path.exists(os.path.join(folder, 'analysis')):
os.mkdir(os.path.join(folder, 'analysis'))
np.save(os.path.join(folder, 'analysis', 'frame_times.npy'), frame_times)
########################################################################################
fig, ax = plt.subplots()
ax.plot(led_vals)
ax.plot(led_frame, np.ones(len(led_frame))*100, '.', color='firebrick')
########################################################################################
fig, ax = plt.subplots()
ax.plot(led_idx / sr, np.ones(len(led_idx)), '.', color='k')
ax.plot(led_frame_to_idx / sr, np.ones(len(led_frame_to_idx))+.1, '.', color='firebrick')
ax.plot([times[0], times[0]], [0.5, 1.5], 'k', lw=1)
ax.plot([times[-1], times[-1]], [0.5, 1.5], 'k', lw=1)
ax.plot(frame_times, np.ones(len(frame_times))*0.5)
ax.set_ylim(0, 2)
plt.show()
embed()
quit()
pass
if __name__ == '__main__':
main(sys.argv[1])

110
event_videos.py Normal file
View File

@ -0,0 +1,110 @@
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import os
import sys
import cv2
import glob
import argparse
from IPython import embed
from tqdm import tqdm
def main(folder, dt):
video_path = glob.glob(os.path.join(folder, '2022*.mp4'))[0]
create_video_path = os.path.join(folder, 'rise_video')
if not os.path.exists(create_video_path):
os.mkdir(create_video_path)
video = cv2.VideoCapture(video_path) # was 'cap'
fish_freqs = np.load(os.path.join(folder, 'analysis', 'fish_freq_interp.npy'))
max_freq, min_freq = np.nanmax(fish_freqs), np.nanmin(fish_freqs)
rise_idx = np.load(os.path.join(folder, 'analysis', 'rise_idx.npy'))
frame_times = np.load(os.path.join(folder, 'analysis', 'frame_times.npy'))
times = np.load(os.path.join(folder, 'times.npy'))
#######################################
for fish_nr in np.arange(2)[::-1]:
for idx_oi in tqdm(np.array(rise_idx[fish_nr][~np.isnan(rise_idx[fish_nr])], dtype=int)):
# idx_oi = int(rise_idx[1][10])
time_oi = times[idx_oi]
# embed()
# quit()
HH = int((time_oi / 3600) // 1)
MM = int((time_oi - HH * 3600) // 60)
SS = int(time_oi - HH * 3600 - MM * 60)
frames_oi = np.arange(len(frame_times))[np.abs(frame_times - time_oi) <= dt]
idxs_oi = np.arange(len(times))[np.abs(times - time_oi) <= dt*3]
fig = plt.figure(figsize=(20/2.54, 20/2.54))
gs = gridspec.GridSpec(2, 1, left=0.1, bottom = 0.1, right=0.95, top=0.95, height_ratios=(4, 1))
ax = []
ax.append(fig.add_subplot(gs[0, 0]))
ax.append(fig.add_subplot(gs[1, 0]))
ax[1].plot(times[idxs_oi] - time_oi, fish_freqs[0][idxs_oi], marker='.', color='firebrick')
ax[1].plot(times[idxs_oi] - time_oi, fish_freqs[1][idxs_oi], marker='.', color='cornflowerblue')
ax[1].set_ylim(min_freq - (max_freq-min_freq)*0.25, max_freq + (max_freq-min_freq)*0.25)
ax[1].set_xlim(-dt*3, dt*3)
ax[0].set_xticks([])
ax[0].set_yticks([])
ax[1].tick_params(labelsize=12)
ax[1].set_xlabel('time [s]', fontsize=14)
# plt.ion()
for i in tqdm(np.arange(len(frames_oi))):
video.set(cv2.CAP_PROP_POS_FRAMES, int(frames_oi[i]))
ret, frame = video.read()
if i == 0:
img = ax[0].imshow(frame)
line, = ax[1].plot([frame_times[frames_oi[i]] - time_oi, frame_times[frames_oi[i]] - time_oi],
[min_freq - (max_freq-min_freq)*0.25, max_freq + (max_freq-min_freq)*0.25],
color='k', lw=1)
else:
img.set_data(frame)
line.set_data([frame_times[frames_oi[i]] - time_oi, frame_times[frames_oi[i]] - time_oi],
[min_freq - (max_freq-min_freq)*0.25, max_freq + (max_freq-min_freq)*0.25])
# label = ('rise_video/frame%4.f.jpg' % len(glob.glob('rise_video/*.jpg'))).replace(' ', '0')
label = (os.path.join(create_video_path, 'frame%4.f.jpg' % len(glob.glob(os.path.join(create_video_path, '*.jpg'))))).replace(' ', '0')
plt.savefig(label, dpi=300)
# plt.pause(0.001)
win_lose_str = 'lose' if fish_nr == 1 else 'win'
# video_name = ("./rise_video/%s_%2.f:%2.f:%2.f.mp4" % (win_lose_str, HH, MM, SS)).replace(' ', '0')
# command = "ffmpeg -r 25 -i './rise_video/frame%4d.jpg' -vf 'pad=ceil(iw/2)*2:ceil(ih/2)*2' -vcodec libx264 -y -an"
video_name = os.path.join(create_video_path, ("%s_%2.f:%2.f:%2.f.mp4" % (win_lose_str, HH, MM, SS)).replace(' ', '0'))
command1 = "ffmpeg -r 25 -i"
frames_path = '"%s"' % os.path.join(create_video_path, "frame%4d.jpg")
command2 = "-vf 'pad=ceil(iw/2)*2:ceil(ih/2)*2' -vcodec libx264 -y -an"
os.system(' '.join([command1, frames_path, command2, video_name]))
os.system(' '.join(['rm', os.path.join(create_video_path, '*.jpg')]))
# os.system(' '.join([command, video_name]))
# os.system('rm ./rise_video/*.jpg')
plt.close()
embed()
quit()
###############################
fig, ax = plt.subplots()
for i, c in enumerate(['firebrick', 'cornflowerblue']):
ax.plot(times, fish_freqs[i], marker='.', color=c)
r_idx = np.array(rise_idx[i][~np.isnan(rise_idx[i])], dtype=int)
ax.plot(times[r_idx], fish_freqs[i][r_idx], 'o', color='k')
pass
##############################
if __name__ == '__main__':
parser = argparse.ArgumentParser(description='Generate videos around events.')
parser.add_argument('file', type=str, help='folder/dataset to generate videos from.')
parser.add_argument('-t', type=float, default=10, help='video duration before and after event.')
# parser.add_argument("-c", action="store_true", help="check if LED pos is correct")
# parser.add_argument('-x', type=int, nargs=2, default=[1272, 1282], help='x-borders of LED detect area (in pixels)')
# parser.add_argument('-y', type=int, nargs=2, default=[1500, 1516], help='y-borders of LED area (in pixels)')
args = parser.parse_args()
main(args.file, args.t)

260
trail_analysis.py Normal file
View File

@ -0,0 +1,260 @@
import os
import sys
import argparse
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import numpy as np
import pandas as pd
from thunderfish.eventdetection import detect_peaks
from IPython import embed
class Trial(object):
def __init__(self, folder, base_path, meta, fish_count):
self._isValid = False
self.base_path = base_path
self.folder = folder
self.meta = meta
self.fish_count = fish_count
self.light_sec = 3 * 60 * 60
self.ids = None
self.fish_freq = None
self.fish_freq_interp = None
self.fish_freq_val = None
self.baseline_freq_times = None
self.baseline_freqs = None
self.rise_idxs = []
self.rise_size = []
self.fish_sign = None
self.fish_sign_interp = None
self.winner = None
self.loser = None
if os.path.exists(os.path.join(self.base_path, self.folder, 'fund_v.npy')):
self.load()
def __repr__(self):
return f'Trial(Date={self.folder}, winner={self.winner})'
# return self.folder
def load(self):
self.fund_v = np.load(os.path.join(self.base_path, self.folder, 'fund_v.npy'))
self.idx_v = np.load(os.path.join(self.base_path, self.folder, 'idx_v.npy'))
self.times = np.load(os.path.join(self.base_path, self.folder, 'times.npy'))
self.ident_v = np.load(os.path.join(self.base_path, self.folder, 'ident_v.npy'))
self.sign_v = np.load(os.path.join(self.base_path, self.folder, 'sign_v.npy'))
self.ids = np.unique(self.ident_v[~np.isnan(self.ident_v)])
if len(self.ids) == self.fish_count:
self.isValid = True
def reshape_and_interpolate(self):
self.fish_freq = np.full((self.fish_count, len(self.times)), np.nan)
self.fish_sign = np.full((self.fish_count, len(self.times), self.sign_v.shape[1]), np.nan)
for enu, id in enumerate(self.ids):
self.fish_freq[enu][self.idx_v[self.ident_v == id]] = self.fund_v[self.ident_v == id]
self.fish_sign[enu][self.idx_v[self.ident_v == id]] = self.sign_v[self.ident_v == id]
self.fish_freq_interp = np.full(self.fish_freq.shape, np.nan)
self.fish_sign_interp = np.full(self.fish_sign.shape, np.nan)
for enu, id in enumerate(self.ids):
i0, i1 = self.idx_v[self.ident_v == id][0], self.idx_v[self.ident_v == id][-1]
self.fish_freq_interp[enu, i0:i1+1] = np.interp(self.times[i0:i1+1],
self.times[self.idx_v[self.ident_v == id]],
self.fish_freq[enu][~np.isnan(self.fish_freq[enu])])
help_sign_v = list(map(lambda x: np.interp(self.times[i0:i1+1], self.times[self.idx_v[self.ident_v == id]], x), self.fish_sign[enu][~np.isnan(self.fish_freq[enu])].T))
self.fish_sign_interp[enu, i0:i1+1] = np.array(help_sign_v).T
def baseline_freq(self, bw = 300):
bins = np.arange(-bw / 2, self.times[-1] + bw / 2, bw)
self.baseline_freq_times = np.array(bins[:-1] + (bins[1] - bins[0])/2)
self.baseline_freqs = np.full((2, len(self.baseline_freq_times)), np.nan)
for enu, id in enumerate(self.ids):
for i in range(len(bins) - 1):
Cf = self.fish_freq[enu][(self.times > bins[i]) & (self.times <= bins[i + 1])]
if len(Cf) == 0:
continue
else:
self.baseline_freqs[enu][i] = np.nanpercentile(Cf, 5)
self.fish_freq_val = [np.nanmean(x[self.baseline_freq_times > self.light_sec]) for x in self.baseline_freqs]
def winner_detection(self):
day_mask = self.times > self.light_sec
day_idxs = np.arange(len(self.times))[day_mask]
shelter_power = np.empty((2, len(day_idxs)))
for enu, id in enumerate(self.ids):
shelter_power[enu] = self.fish_sign_interp[enu][day_idxs, -1]
mean_shelter_power = np.nanmean(shelter_power, axis=1)
self.winner = 1 if mean_shelter_power[1] > mean_shelter_power[0] else 0
self.loser = 0 if self.winner == 1 else 1
def rise_detection(self, rise_th):
def check_rises_size(peak):
peak_f = self.fish_freq[i][peak]
peak_t = self.times[peak]
closest_baseline_idx = list(map(lambda x: np.argmin(np.abs(self.baseline_freq_times - x)), peak_t))
closest_baseline_freq = self.baseline_freqs[i][closest_baseline_idx]
rise_size = peak_f - closest_baseline_freq
return rise_size
for i in range(len(self.fish_freq)):
rise_peak_idx, trough = detect_peaks(self.fish_freq[i][~np.isnan(self.fish_freq[i])], rise_th)
non_nan_idx = np.arange(len(self.fish_freq[i]))[~np.isnan(self.fish_freq[i])]
rise_peak_idx, trough = non_nan_idx[rise_peak_idx], non_nan_idx[trough]
rise_size = check_rises_size(rise_peak_idx)
self.rise_idxs.append(rise_peak_idx[rise_size >= rise_th])
self.rise_size.append(rise_size[rise_size >= rise_th])
def update_meta(self):
entries = self.meta.index.tolist()
if self. folder not in entries:
self.meta.loc[self.folder] = ['' for _ in self.meta.columns]
self.meta.loc[self.folder, 'Win_ID'] = self.ids[self.winner]
self.meta.loc[self.folder, 'Lose_ID'] = self.ids[self.loser]
self.meta.loc[self.folder, 'Win_EODf'] = self.fish_freq_val[self.winner]
self.meta.loc[self.folder, 'Lose_EODf'] = self.fish_freq_val[self.loser]
self.meta.loc[self.folder, 'Win_rise_c'] = len(self.rise_idxs[self.winner])
self.meta.loc[self.folder, 'Lose_rise_c'] = len(self.rise_idxs[self.loser])
self.meta.loc[self.folder, 'light_sec'] = self.light_sec
self.meta.to_csv(os.path.join(self.base_path, 'meta.csv'), sep =',')
def ilustrate(self):
fig = plt.figure(figsize=(20/2.54, 12/2.54))
gs = gridspec.GridSpec(1, 1, left = 0.1, bottom = 0.1, right = 0.95, top = 0.95)
ax = fig.add_subplot(gs[0, 0])
for enu, id in enumerate(self.ids):
c = 'firebrick' if self.winner == enu else 'forestgreen'
ax.plot(self.times, self.fish_freq[enu], marker='.', color=c, zorder=1)
ax.plot(self.times[np.isnan(self.fish_freq[enu])], self.fish_freq_interp[enu][np.isnan(self.fish_freq[enu])], '.', zorder=1, color=c, alpha=0.25)
ax.plot(self.baseline_freq_times, self.baseline_freqs[enu], '--', color='k', zorder=2)
ax.plot(self.times[self.rise_idxs[enu]], self.fish_freq_interp[enu][self.rise_idxs[enu]], 'o', color='k')
win_str = '(W)' if self.winner == enu else ''
ax.text(self.times[-1], self.fish_freq_val[enu]-10, '%.0f' % id + win_str, va ='center', ha='right')
ax.set_xlim(0, self.times[-1])
freq_range = (np.nanmin(self.fish_freq), np.nanmax(self.fish_freq))
ax.set_ylim(freq_range[0] - 20, freq_range[1] + 10)
plt.show()
def save(self):
saveorder = -1 if self.winner == 1 else 1
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'fish_freq.npy'), self.fish_freq[::saveorder])
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'fish_freq_interp.npy'), self.fish_freq_interp[::saveorder])
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'baseline_freqs.npy'), self.baseline_freqs[::saveorder])
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'baseline_freq_times.npy'), self.baseline_freq_times[::saveorder])
help_lens = [len(x) for x in self.rise_idxs]
rise_idxs_s = np.full((self.fish_count, np.max(help_lens)), np.nan)
rise_size_s = np.full((self.fish_count, np.max(help_lens)), np.nan)
for i in range(self.fish_count):
rise_idxs_s[i][:len(self.rise_idxs[i])] = self.rise_idxs[i]
rise_size_s[i][:len(self.rise_size[i])] = self.rise_size[i]
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'rise_idx.npy'), rise_idxs_s[::saveorder])
np.save(os.path.join(self.base_path, self.folder, 'analysis', 'rise_size.npy'), rise_size_s[::saveorder])
@property
def isValid(self):
return self._isValid
@isValid.setter
def isValid(self, value):
print('Trial (%s) is valid' % (self.folder))
self._isValid = value
def frame_to_idx(self, event_frames):
self.sr = 20000
LED_idx = pd.read_csv(os.path.join(self.folder, 'led_idxs.csv'), sep=',')
led_idx = np.array(LED_idx).T[0]
led_frame = np.load(os.path.join(self.folder, 'LED_frames.npy'))
led_idx_span = led_idx[-1] - led_idx[0]
led_frame_span = led_frame[-1] - led_frame[0]
frames_to_idx = ((event_frames - led_frame[0]) / led_frame_span) * led_idx_span + led_idx[0]
event_times = frames_to_idx / self.sr
return event_times
def main():
parser = argparse.ArgumentParser(description='Evaluated electrode array recordings with multiple fish.')
parser.add_argument('-f', type=str, help='single recording analysis', default='')
# parser.add_argument("-c", action="store_true", help="check if LED pos is correct")
# parser.add_argument('-x', type=int, nargs=2, default=[1272, 1282], help='x-borders of LED detect area (in pixels)')
# parser.add_argument('-y', type=int, nargs=2, default=[1500, 1516], help='y-borders of LED area (in pixels)')
args = parser.parse_args()
base_path = '/home/raab/data/2022_competition'
if os.path.exists(os.path.join(base_path, 'meta.csv')):
meta = pd.read_csv(os.path.join(base_path, 'meta.csv'), sep=',', index_col=0)
else:
meta = None
if args.f == '':
folders = os.listdir(base_path)
folders = [x for x in folders if not '.' in x]
else:
folders= [os.path.split(os.path.normpath(args.f))[-1]]
trials = []
for folder in folders:
trial = Trial(folder, base_path, meta, fish_count=2)
if not trial.isValid:
continue
trial.reshape_and_interpolate()
trial.winner_detection()
trial.baseline_freq(bw=300)
# ToDo: q10 corrected EODfs
trial.rise_detection(rise_th=5)
if meta is not None:
trial.update_meta()
trial.save()
trial.ilustrate()
trials.append(trial)
# meta.loc[folder, 'Fish1_ID'] = 1
# meta.to_csv('')
if __name__ == '__main__':
main()