Skip to content
Snippets Groups Projects
Select Git revision
  • 755203c149a42e95449da95b1e516337126a7ac5
  • main default protected
2 results

activity.jl

Blame
  • activity.jl 3.50 KiB
    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.2)
    	frequency_activity = 1*(frequency.>2.2)+1*(frequency.>1.0)
    	
    
    	# Get aggregate activity:
    	activity = zeros(length(power_activity),1) 
    	activity[1] = power_activity[1] #magnitude data looks slightly more reliable, especially when standing still
    	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