diff --git a/src/fox_lab_to_nwb/assets/trialanalysis_photron.m b/src/fox_lab_to_nwb/assets/trialanalysis_photron.m
new file mode 100644
index 0000000..d2ba1aa
--- /dev/null
+++ b/src/fox_lab_to_nwb/assets/trialanalysis_photron.m
@@ -0,0 +1,383 @@
+%INPUT: single trial folder, plot (0 or 1) 
+%analyzes all data within a file and saves struct as .pcf
+%can use this per trial, called as part of "folderanalysis"
+%toplot and tosave should =1 if yes
+
+function trialdata = trialanalysis_photron(trial_path, env, toplot, tosave)
+
+%if not all the arguments, redefine
+%we ain't about checking for special cases
+if nargin < 3
+    trial_path = uigetdir();
+    toplot = 1;
+    env = 40; %CHANGE ENVELOPE SIZE HERE if running single trial
+    tosave = 0;
+end
+
+trialdata = struct;
+
+%load in data
+%daq
+fd = dir([trial_path filesep '*.fly2']);
+ddata = load(fullfile(fd.folder, filesep, fd.name),'-mat');
+ddata = ddata.rec;
+
+%haltere data
+hd = dir([trial_path filesep '*_haltamp.mat']);
+hdata = load(fullfile(hd.folder, filesep, hd.name));
+
+%antenna/head data
+ad = dir([trial_path filesep '*_antennatracking*.csv']);
+%make a variable to know if head fixed or head free
+% headfree = 0;
+% if isempty(ad)
+%     ad = dir([trial_path filesep '*_headtracking*.csv']);
+%     headfree = 1;
+% end
+adata = readtable(fullfile(ad.folder, filesep, ad.name), 'NumHeaderLines',1);
+
+%wing data
+wd = dir([trial_path filesep '*_PROC.mat']);
+wdata = load(fullfile(wd.folder, filesep, wd.name));
+
+%antenna columns: 
+%1 = coords (frame number), then x y likelihood for each
+%2-4 = L_ant; 5-7 = R_ant; 8-10 = L_base; 11-13 = R_base
+%head tracking:
+%14-16 = L_OC; 17-19 = R_OC; 20-22 = L_pOR; 23-25 = R_pOR; 
+%26-28 = L_aVT; 29-31 = R_aVT
+
+%ASSIGN DATA VARIABLES
+%daq
+%trial info
+duration = ddata.trial.duration;
+stimlength = ddata.trial.stimlength;
+flynum = ddata.trial.flynum;
+trialnum = ddata.trial.flytrial;
+volt = ddata.trial.stimvolt;
+%data
+stamps = ddata.daq.tstamps;
+ctrig = ddata.daq.data(:,2); %fastec cam trigger
+%otrig = ddata.daq.data(:,3); %opto trigger
+ptrig = ddata.daq.data(:,9); %photron trigger
+wtrig = ddata.daq.data(:,10); %wind trigger
+
+%NOT CURRENTLY USING
+%including in struct though so it's there if needed?
+%otrig = ddata.daq.data(:,3); %opto trigger
+dwbf = ddata.daq.data(:,6);
+%wbf data from wingbeat amplifier needs to be multiplied by 100: 
+% 1V = 100Hz
+dwbf = dwbf*100; %daq wingbeat freq, added 7/17
+dwbaL = ddata.daq.data(:,4); %temp
+% twbaL = twbaL - mean(twbaL(1:2499)); %set to 0 for denoising
+dwbaR = ddata.daq.data(:,5);
+% twbaR = twbaR - mean(twbaR(1:2499));
+% hutchenL = ddata.daq.data(:,7); 
+% hutchenR = ddata.daq.data(:,8); 
+
+%ALIGN CAMERAS TO DAQ - CHECK WITH MIKE AGAIN?
+%find trigger times
+ctrigtime = find(ctrig>3,1)/ddata.daq.fs;
+ptrigtime = find(ptrig>3,1)/ddata.daq.fs;
+
+%topcam (currently skipping sidecam!)
+meta_t = fastecMetaReader(fullfile(trial_path, 'TOPCAM_000000.txt'));
+ts_t = linspace(1/meta_t.fs, meta_t.numframes/meta_t.fs, meta_t.numframes);
+ts_t = ts_t-(ts_t(end)-ctrigtime);
+%photron/phantom camera
+mp = dir([trial_path filesep 'HALTCAM*.mii']); %since the number changes every time
+meta_p = photronMetaReader(fullfile(mp.folder, mp.name));
+ts_p = linspace(1/meta_p.fs, meta_p.numframes/meta_p.fs, meta_p.numframes);
+ts_p = ts_p-(ts_p(end)-ptrigtime);
+
+%% ANTENNAE
+%subtract height from y data because y = 0 is at top of image
+Lx = adata.L_ant;
+Ly = meta_t.height - adata.L_ant_1;
+Rx = adata.R_ant;
+Ry = meta_t.height - adata.R_ant_1;
+%base data just in case
+Lxb = adata.L_base;
+Lyb = meta_t.height - adata.L_base_1;
+Rxb = adata.R_base;
+Ryb = meta_t.height - adata.R_base_1;
+%find theta
+Lt = atan2d((Ly-Lyb),(Lx-Lxb));
+Lt = wrapTo360(Lt);
+Rt = atan2d((Ry-Ryb),(Rx-Rxb));
+%find distance between L and R
+dant = [];
+for j = 1:length(Lx)
+    dant(j) = pdist2([Lx(j) Ly(j)],[Rx(j) Ry(j)]);
+end  
+dant = dant';
+
+%% WING DATA
+%denoise daq WBA freq data
+dwbf = hampel(dwbf);
+%flyalyzer data
+wbL = wdata.fly.wingL.angle;
+% wbLr = wdata.fly.wingL.root;
+wbR = wdata.fly.wingR.angle;
+%need to wrap R wing data
+wbR = wrapTo180(wbR);
+% wbRr = wdata.fly.wingR.root;
+%create timestamps: flyalyzer tstamps aren't good, it's always every 20 
+wstamps = linspace(ts_t(1),ts_t(end),length(ts_t)/20);
+
+%zeroing removed 7/17: will be zeroed later!
+%zero everything, camera upix is 1351 -> 1350/20 = 67.5 - see
+%"cameratimes_phantom.m" for calculation
+% wbL = wbL-mean(wbL(1:60));
+% wbR = wbR-mean(wbR(1:60));
+
+%7/17 alright we're gonna modify the freq calc bc something's wrong
+%NOTE: Mike's code, this sums up image pixels to get the frequency. This
+%does include the haltere because it's in the video, but the wings are much
+%bigger so it shouldn't make too much difference.
+%MODIFIED SUMVID
+vr = VideoReader([trial_path filesep 'TOPCAM_000000.avi']);
+%create array for video sums
+int = zeros(1,vr.NumFrames);
+%go through frames
+for j = 1:vr.NumFrames
+    %read image, greyscale but saves as color
+    img = read(vr,j);
+    %take only red channel since they're all the same
+    img = img(:,:,1);
+    %add up all the pixels in that image and add to time series
+    int(j)=sum(sum(img));
+end
+%pxx is made of FFTs for each timestamp, f is freq, t is time
+%zscore gets rid of DC info (low freq)
+%window size is 64 w 63 overlap (or 128/127), or use 20/0 so same
+%sample number as flyalyzer
+%freq range is 150-300, 2000 is fps, 'reassigned' does a thing
+%reassigned is kinda like contrast stretching?
+[pxx,f,t]=spectrogram(zscore(int),64,63,150:300,2000,'reassigned');
+%time/freq ridge is pretty much the wbf plot, has a .001 penalty to
+%make sure there's not much difference between y values (so it doesn't
+%jump)
+%NOTE: there might be some weird harmonics if the penalty isn't small
+%enough, can check by plotting t,calc_wbf 
+calc_wbf = tfridge(pxx,f,.001);
+
+%NOTE: there might be some weird harmonics if the penalty isn't small
+%enough, can check by plotting t,calc_wbf over top of the spectrogram to
+%check if it lines up
+    % spectrogram(zscore(int),128,127,150:300,2000,'yaxis','reassigned');
+    % plot(t,calc_wbf,'r')
+%add NaNs back in to beginning so same length as timestamps
+calc_wbf = [NaN((length(int) - length(calc_wbf)),1); calc_wbf];
+
+%% Haltere Data
+%in ampalyzer, haltere base needs to be to the right of the haltere
+%base, not in the center! otherwise lots of weird artefacts
+%haltere data
+%get coordinates and subtract y from height of frame (y=0 at top)
+p_height = meta_p.height;
+%haltere and wing root
+hroot = hdata.haltrpos;
+hroot(2) = p_height-hroot(2);
+wroot = hdata.wingrpos;
+wroot(2) = p_height-wroot(2);   
+hx = hdata.haltpos(:,1);
+hy = p_height - hdata.haltpos(:,2);
+%this bit is from Mike and some trial and error...
+%find a transformation matrix that turns xy in haltere root to 0,0 and xy in
+%wing root to 1,0
+tform = fitgeotrans([hroot;wroot],[0,0;1,0],'nonreflectivesimilarity');
+%use that transformation matrix for the tracked points
+[hxn, hyn] = transformPointsForward(tform,hx,hy);
+%get angles and wrap
+Ht = atan2d(hyn,hxn); %PHOTRON
+%get AMPLITUDE (envelope function)
+%envelope of 14 looks best for phantom, about 40 for photron -DEFINED AT TOP THIS IS WHY
+%WE'RE GOING THROUGH FOLDERS DAMMIT
+[Htup,Htlow] = envelope(Ht,env,'peak');
+%find peaks the other way... 10k fps at 200 Hz peak-to-peak should be 50,
+%so half a stroke is 25, should be just above that...
+[Htu_pks, Htu_locs, Htu_w] = findpeaks(Ht, 'MinPeakHeight',20, 'MinPeakDistance',30);
+[Htl_pks, Htl_locs, Htl_w] = findpeaks(-Ht, 'MinPeakHeight',20, 'MinPeakDistance',30);
+Htl_pks = -Htl_pks;
+% %distance, normalized to distance between wing and root
+% %removed 6/28 bc not a great time? kept Hd though for reference
+% %alternates between top and bottom of wing stroke; every other peak has a
+% %hutchen
+Hd = pdist2([0 0],[hxn,hyn]);
+% %dips at top of stroke....
+% [Hdup,Hdlow] = envelope(Hd,env,'peak');
+% [Hdu_pks, Hdu_locs, Hdu_w] = findpeaks(Hd, 'MinPeakHeight',.4, 'MinPeakDistance',30);
+
+%get FREQUENCY
+%spectrogram function basically makes an fft (power x function: pxx) for
+%each time t, samples every 256 points (255 overlapping)
+% 64 & 128 is too small because fps is 10k! 
+%to be optimal range
+%detrend gets rid of all the low freq info (dc info?)
+%looking for first set of freq and no harmonics (max wbf is around 230?)
+[pxx,f,~]=spectrogram(detrend(Ht),256,255,0:300,meta_p.fs);
+%find time frequency ridge (kinda like the avg?)
+Hfreq = tfridge(pxx,f,.001);
+%since spec uses sliding window, first set of overlap is not used, add NaNs
+%to make same size as timestamps
+Hfreq = [NaN((length(Ht) - length(Hfreq)),1); Hfreq];
+
+if isequal(toplot,1) %plot data
+    close all
+    figure('WindowState','maximized');
+    t =  tiledlayout(2,3);
+
+    nexttile
+    hold on
+    %yyaxis left
+    plot(wstamps,wbR,'r')
+    %for some reason this is dashed if you don't specify it...
+    plot(wstamps,wbL,'-b') 
+    % yyaxis right
+    % plot(stamps,wbf)
+    title('Wingbeat Amp R (red) & L (blue)')
+    hold off
+    
+    %7/17 look at freq
+    nexttile
+    hold on
+    yyaxis left
+    ylabel('WB Amp')
+    plot(stamps,dwbf)
+    yyaxis right
+    ylabel('Calc')
+    plot(ts_t,calc_wbf)
+    
+    hold off
+    nexttile
+    hold on
+    plot(ts_t,dant)
+    title('Distance btwn Antennae')
+    hold off
+
+    nexttile
+    hold on
+    yyaxis left
+    plot(ts_t,Lt,'b')
+    ylabel("Ltheta")
+    yyaxis right
+    plot(ts_t,Rt,'r')
+    ylabel("Rtheta")
+    title('Antennal angle')
+    hold off
+
+    nexttile
+    hold on
+    plot(ts_p,Hfreq)
+    title('Haltere Freq')
+    hold off
+
+    % nexttile
+    % hold on
+    % yyaxis left
+    % plot(ts_p,(Htup-Htlow),'r')
+    % ylim([150 200])
+    % yyaxis right
+    % plot(ts_p,(Hdup-Hdlow),'b')
+    % ylim([.35 .65])
+    % title('Haltere Envelope Change angle (red) & length (blue)')
+    % hold off
+
+    nexttile
+    hold on
+    envelope(Ht,env,'peak')
+    title('Haltere angle')
+    hold off
+
+    % nexttile
+    % hold on
+    % plot(Ht)
+    % plot(Htu_locs,Htu_pks,'o')
+    % plot(Htl_locs,Htl_pks,'o')
+    % title('Haltere angle peaks')
+    % hold off
+
+
+    % nexttile
+    % hold on
+    % envelope(Hd,env,'peak')
+    % plot(Hdu_locs,Hdu_pks,'o')
+    % plot(Hdl_locs,Hdl_pks,'o')
+    % title('Haltere dist')
+    % hold off
+
+    title(t,{"Fly: " + num2str(flynum) + " Trial: " + num2str(trialnum), ...
+        "Stim " + num2str(stimlength) + " Voltage: " + num2str(volt), ...
+        "Envelope: " + num2str(env)});
+end
+
+ %% Add to struct
+ %note: nothing is normalized to the mean yet! This is only the raw
+ %data!
+ %not added: hutchens
+ %trial info
+ trialdata.info.duration = duration;
+ trialdata.info.stimlength = stimlength;
+ trialdata.info.flynum = flynum;
+ trialdata.info.trialnum = trialnum;
+ trialdata.info.envelope = env;
+ trialdata.info.volt = volt;
+ %daq data
+ trialdata.trial.tstamps = stamps;
+ trialdata.trial.ctrig = ctrig;
+ %trialdata.trial.otrig = otrig;
+ trialdata.trial.wbaL = dwbaL;
+ trialdata.trial.wbaR = dwbaR;
+ trialdata.trial.wbf = dwbf;
+ trialdata.trial.ptrig = ptrig;
+ trialdata.trial.wtrig = wtrig;
+ %camera data
+ trialdata.ant.tstamps = ts_t';
+ trialdata.ant.dist_ant = dant;
+ trialdata.ant.Rt = Rt;
+ trialdata.ant.Lt = Lt;
+ trialdata.wing.wbL = wbL;
+ trialdata.wing.wbR = wbR;
+ trialdata.wing.wstamps = wstamps;
+ trialdata.wing.wbf = calc_wbf;
+ trialdata.halt.tstamps = ts_p';
+ trialdata.halt.Ht = Ht;
+ trialdata.halt.Htup = Htup;
+ trialdata.halt.Htlow = Htlow;
+ trialdata.halt.Hdist = Hd;
+ % trialdata.halt.Hdup = Hdup;
+ % trialdata.halt.Hdlow = Hdlow;
+ trialdata.halt.Hfreq = Hfreq;
+ trialdata.halt.Htu_pks = Htu_pks;
+ trialdata.halt.Htu_locs = Htu_locs;
+ trialdata.halt.Htu_w = Htu_w;
+ trialdata.halt.Htl_pks = Htl_pks;
+ trialdata.halt.Htl_locs = Htl_locs;
+ trialdata.halt.Htl_w = Htl_w;
+ % trialdata.halt.Hdu_pks = Hdu_pks;
+ % trialdata.halt.Hdu_locs = Hdu_locs;
+ % trialdata.halt.Hdu_w = Hdu_w;
+ % trialdata.halt.Hdl_pks = Hdl_pks;
+ % trialdata.halt.Hdl_locs = Hdl_locs;
+ % trialdata.halt.Hdl_w = Hdl_w;
+
+ % %head data if applicable TO DO
+ % if headfree == 1
+ %     trialdata.head.pitch = [];
+ %     trialdata.head.yaw = [];
+ %     trialdata.head.roll = [];
+ % end
+
+ if tosave == 1
+     %% save
+     %save analysis so it runs faster if you're playing around with the graphs
+     %date & time
+     %dt = string(datetime('now'),'yyMMdd_HHmmss');
+     %save file name (Tsh file)
+     flycross = extractBetween(trial_path,'trials\','\');
+     flycross = convertCharsToStrings(flycross);
+     sfn = convertCharsToStrings(trial_path) + filesep + flycross + '_f' + flynum + '_r' + trialnum + '_trial.pcf';
+     save(sfn,'trialdata');
+ end
\ No newline at end of file
diff --git a/src/fox_lab_to_nwb/assets/windtrial.m b/src/fox_lab_to_nwb/assets/windtrial.m
new file mode 100644
index 0000000..de9bcb5
--- /dev/null
+++ b/src/fox_lab_to_nwb/assets/windtrial.m
@@ -0,0 +1,181 @@
+function rec = windtrial
+%AMS - edited from "windtrial_phantom" in "Antennae_Opto_0623"
+%daq setup is for opto trials - may need to modify later
+%to use with new photron camera
+
+%returns struct (rec) as fly2 file
+
+%need to reset daq, don't need to clear all b/c you're in a function
+daqreset
+
+%% DEFINE VARIABLES HERE FOR EACH TRIAL
+%fly info for save file   
+flynum = 16;
+flytrial = 1;    
+flycross = 'TESTING'; 
+%directory name UPDATE IF NEEDED
+datafolder = "C:\Users\Foxlab\Dropbox\Fox_lab\Amy\Antennae_windtrials";
+
+%stimulus info (in seconds) 
+
+stimtime = .5; %length of stimulus 
+stimint = 1; %interval between stimuli/after stimulus
+stimwait = 0.25; %time before first stim
+stimrep = 1; %repetitions 
+
+%leaving opto stuff in FOR NOW because of how daq is set up!
+%set stimvolt to 0 so it doesn't go off
+stimvolt = 0; 
+
+%% END VARIABLE EDITS
+
+
+%% MESSAGE BOX 
+%this is just copy/pasted from the help page for the inputdlg function
+%not as good as a gui like xdrift b/c doesn't save overwritten variables
+prompt = {'Cross', 'Fly Number','Trial Number', 'Stimulus Length', 'Inter-Stimulus Interval',...
+    'Lag Time', 'Repetitions', 'Opto Volt'};
+dlgtitle = 'Trial Data';
+fieldsize = [1 30; 1 30; 1 30; 1 30; 1 30; 1 30; 1 30; 1 30];
+%default values
+definput = {flycross, num2str(flynum), num2str(flytrial), num2str(stimtime), num2str(stimint),...
+    num2str(stimwait), num2str(stimrep), num2str(stimvolt)};
+%creates a cell array
+answer = inputdlg(prompt,dlgtitle,fieldsize,definput);
+%if there's something empty then cancel trial
+if isempty(answer)
+    return
+end
+
+%rewrite data
+flycross = answer{1};
+flynum = str2double(answer{2});
+flytrial = str2double(answer{3});
+stimtime = str2double(answer{4});
+stimint = str2double(answer{5});
+stimwait = str2double(answer{6});
+stimrep = str2double(answer{7});
+stimvolt = str2double(answer{8});
+
+
+%define total stim length
+stimlength = stimwait + stimrep*(stimtime+stimint); 
+
+%% Make DAQ objects & channels  
+%INITIATE DAQ
+d = daq('ni'); 
+d.Rate = 10000; %10k rate
+
+%ADD DAQ CHANNELS
+fastectrig = addoutput(d, "Dev2","ao0","Voltage");
+optotrig = addoutput(d,"Dev2","ao1","Voltage"); %DISCONNECTED
+phantomtrig = addoutput(d,"Dev2","ao2","Voltage");
+windtrig = addoutput(d, "Dev2","ao3","Voltage"); %WIND
+
+CamSync = addinput(d, "Dev2", "ai0","Voltage"); 
+CamTrigger = addinput(d, "Dev2", "ai1","Voltage"); 
+OptoTrigger = addinput(d, "Dev2", "ai2","Voltage");
+LWingbeatAmp = addinput(d, "Dev2", "ai3","Voltage");
+RWingbeatAmp = addinput(d, "Dev2", "ai4","Voltage");
+WingbeatFreq = addinput(d, "Dev2", "ai5","Voltage");
+LHutchen = addinput(d, "Dev2", "ai6","Voltage");
+RHutchen = addinput(d, "Dev2", "ai7","Voltage");
+%there's no AI 8-15 on the daq?
+PTrigger = addinput(d, "Dev2", "ai16","Voltage");
+WTrigger = addinput(d, "Dev2", "ai17","Voltage");
+
+%START SUB-FUNCTION (this would return "scandata" and not rec!)
+
+%CREATE TRIGGER/INPUT VECTORS
+%camera trigger is the stimlength multiplied by the rate, then UP about
+%.25sec after? could use the stimtime for the trigger but if you want the
+%stimulus to be longer then the file will be too big. 
+%all input vectors need to be single columns, and all the same size!
+%camera trigger - goes up a the end of the trial
+%multiply by 3.3 because the fastec cameras want a voltage of 3.3
+trigtime = 0.25; %time after trigger/when trig is 1! using 250 ms so 
+camtrig = [zeros(stimlength*d.Rate,1); ones(trigtime*d.Rate,1)]*3.3;
+%Phantom camera is TTL! needs 5V to trigger, won't trigger at 3.3
+phantomtrig = [zeros(stimlength*d.Rate,1); ones(trigtime*d.Rate,1)]*5;
+totallength = stimlength + trigtime;
+%KEEPING THIS IN FOR NOW, OPTO LIGHT NOT IN PLACE
+%opto trig, might have variable sizes and lengths so add with a for loop
+otrig = [zeros(stimwait*d.Rate,1)]; %add rate time
+for i=1:stimrep %COULD EDIT
+    otrig = [otrig; ones(stimtime*d.Rate,1); zeros(stimint*d.Rate,1)];
+end
+%need to add to the end so it's the same length as the trigger
+otrig = [otrig; zeros(0.25*d.Rate,1)];
+%multiply by voltage if we want to vary intensity (change LED driver to
+%"mod")
+otrig = otrig*stimvolt;
+
+%WIND TRIGGER
+wtrig = [zeros(stimwait*d.Rate,1)]; %add rate time
+for i=1:stimrep %COULD EDIT
+    wtrig = [wtrig; ones(stimtime*d.Rate,1); zeros(stimint*d.Rate,1)];
+end
+wtrig = [wtrig; zeros(0.25*d.Rate,1)]; 
+windvolt = 9; %solenoid meant to run at 7V, doesn't react until 9V
+%solenoid power source should be 12V btw
+wtrig = wtrig*windvolt; 
+
+
+%indicator LED
+
+%% GO!
+%create output matrix, needs to be MxN matrix, where M is num data scans
+%and N is num output channels
+outdata = [camtrig otrig phantomtrig wtrig];
+%readwrite function: inScanData = readwrite(d,outScanData); 
+scandata = readwrite(d,outdata,"OutputFormat","Matrix");
+
+%% CREATE STRUCT AND SAVE FILE
+%create struct, needs to follow xdrift naming!
+rec = struct;
+%daq data
+rec.daq.fs = d.Rate;
+rec.daq.data = scandata;
+rec.daq.tstamps = linspace(0,totallength,totallength*d.Rate)';
+rec.daq.channelnames = ["CamSync" "CamTrigger" "OptoTrigger" "LWingBeatAmp" "RWingBeatAmp" "WingBeatFreq" "LHutchen" "RHutchen" "PTrigger" "WindTrigger"];
+
+%need to create base file name here for save
+%need to convert flycross from char vector to string
+flycross = convertCharsToStrings(flycross);
+ffn = flycross + '_' + string(datetime('now'),'yyMMdd_HHmmss') + '_f' + string(flynum) + '_r' + string(flytrial);
+
+%trial info, can also be figured out from daq data
+rec.trial.info = ffn; %trial info (filename)
+rec.trial.duration = stimlength;
+rec.trial.stimlength = stimtime;
+rec.trial.interval = stimint;
+rec.trial.delay = stimwait;
+rec.trial.reps = stimrep;
+rec.trial.flynum = flynum;
+rec.trial.flytrial = flytrial;
+rec.trial.stimvolt = stimvolt;
+rec.trial.windvolt = windvolt;
+
+%SAVE
+%create trial folder, this will also create the cross folder if it doesn't
+%exist! just make sure the datafolder is correct!
+%directory file name
+dfn = datafolder + filesep + flycross + filesep + ffn;
+mkdir(dfn);
+%save file name (fly2 file)
+sfn = dfn + filesep + ffn + '.fly2';
+save(sfn,'rec');
+
+%END SUB-FUNCTION HERE?
+
+%PLOT TO CHECK WBA
+wbaL = rec.daq.data(:,4);
+wbaR = rec.daq.data(:,5);
+wbf = rec.daq.data(:,6);
+plot(wbaL - wbaL(1))
+hold on
+plot(wbaR - wbaR(1))
+plot(wbf - wbf(1))
+hold off
+
+end