Finalized version

frequency_powerFINAL.jl

+using Pkg; Pkg.activate(".");
+using FFTW, DSP, Plots, DelimitedFiles
+using REPL.TerminalMenus
+
+# Plot heatmap of spectrogram:
+function plot_spectrogram(s)
+	n = 80 #window size
+	noverlap = n-1 #window overlap
+	spec = spectrogram(s, n, noverlap)
+	power = spec.power./sum(spec.power, dims=1)
+	heatmap(spec.time, spec.freq, power, size=(1500,800),title="Spectrogram")
+end
+
+# Get frequency with greatest power:
+function get_spectrogram(s)
+	n = 80 #window size
+	noverlap = n-1 #window overlap
+	spec = spectrogram(s, n, noverlap)
+	power = spec.power./sum(spec.power, dims=1)
+	power = 0.5*10*maxPower(power)/41 #adjust to get correct unit
+	return power
+end
+
+# Get index of largest value in each column:
+function maxPower(data)
+	l = size(data)[2]
+	index = zeros(l,1)
+	for i=1:l
+		index[i] = findmax(data[:,i])[2]
+	end
+	return index
+end
+
+#Moving mean function for smoothing magnitude of acceleration in time domain:
+movmean(x, n) = conv(x, ones(Int64, n))/n
+
+# Load data:
+function read_data(dataname, sensor) 
+	return readdlm("data/$(dataname)_$(sensor).txt", Float64)
+end
+
+read_acc(dataname) = read_data(dataname, "acc")
+normalize_data(A) = sqrt.(sum(A.^2, dims=2)) #Euclidean norm of data
+decimate(v) = resample(v, 0.1) #decimate data
+mean(v) = sum(v)/length(v)
+detrend(v) = v .- mean(v) #remove mean
+xyz(data) = data[:,2:4] #extract relevant components
+
+# Acceleration data pipieline:
+acc_data(dataname) = read_acc(dataname) |> xyz |> normalize_data |> vec |> decimate
+
+# Spectrogram calculation:
+acc_frequency(dataname) =  acc_data(dataname) |> detrend |> get_spectrogram
+acc_spectrogram(dataname) = acc_data(dataname) |> detrend |> plot_spectrogram
+
+# Moving average of magnitude of acceleration:
+n_smooth = 80 #kernel size
+acc_power(dataname) = acc_data(dataname) |> x->movmean(x,n_smooth) |> x->x[n_smooth:end-n_smooth+1]
+
+function menuplot()
+	datasets = ["jogging1", "jogging2", "mix1", "standing1", "walking1", "walking2"]
+	
+	# Menu:
+	keybindings = [Char('0' + i) for i in 1:length(datasets)]
+	menu = RadioMenu(datasets, pagesize=min(length(datasets), 10), keybindings=keybindings)
+	choice = request("Please choose a dataset (up/down or 1-9):", menu)
+	
+	# Get data:
+	frequency = acc_frequency(datasets[choice]) #get frequency data
+	power = acc_power(datasets[choice]) #get magnitude data
+	plot_heatmap = acc_spectrogram(datasets[choice]) #plot spectrogram heatmap
+
+	# Decide activity:
+	power_activity = 1*(power.>11.5)+1*(power.>10)
+	frequency_activity = 1*(frequency.>2.2)+1*(frequency.>1.0)
+	
+
+	# Get aggregate activity:
+	activity = zeros(length(power_activity),1) 
+	activity[1] = frequency_activity[1] #frequency data looks slightly more reliable
+	for t = 2:length(power_activity)
+		if power_activity[t] == frequency_activity[t]
+			activity[t] = power_activity[t] #update activity if datasets agree
+		else
+			activity[t] = activity[t-1] #otherwise continue previous activity
+		end
+	end
+	
+	# Plotting:
+	plot_PA = plot(power_activity,size=(1500,800),title="Activity (magnitude)")
+	plot_FA = plot(frequency_activity,size=(1500,800),title="Activity (frequency)")	
+	plot_activity = plot(activity,size=(1500,800),title="Activity (aggregate)")
+	plot_frequency = plot(frequency,size=(1500,800),title="Frequency")
+	plot_power = plot(power,size=(1500,800),title="Magnitude")
+	
+	display(plot(plot_heatmap,plot_frequency,plot_power,plot_PA,plot_FA,plot_activity))
+	
+	# Combined plots:
+	#plot(frequency,size=(1500,800),label="Frequency")
+	#plot!(activity,label="Activity")
+	#plot(power.-9.8,size=(1500,800),label="Magnitude")
+	#plot!(activity,label="Activity")
+end
+