Merge 8caa30e into 65582b2

.gitignore

@@ -10,6 +10,12 @@ test/media/
 # media files generated by scripts
 scripts/media/
 
+# ipynb checkpoints in pythonscripts
+scripts/pythonscripts/.ipynb_checkpoints
+
+# idm acc trace files stored in lmfit
+scripts/lmfit/
+
 # Files generated by invoking Julia with --code-coverage
 *.jl.cov
 *.jl.*.cov

scripts/experiments.jl

@@ -86,6 +86,24 @@ rmse_plots_modelscompare(rmse_pos_list,filename = "media/rmse_pos_upper.svg");
 coll_mat_list = [coll_mat_idm,coll_mat_cidm,coll_mat_lmidm,coll_mat_pf];
 coll_barchart(coll_mat_list,filename = "media/coll_barchart_upper.svg");
 
+#*****************Generate idm params for lmfit************
+s = scenarios_upper
+f = FilteringEnvironment()
+name = "upper"
+
+for i in 1:length(s)
+        veh_id_list = s[i][1]
+        ts = s[i][2][1]
+        te = s[i][2][2]
+        filename = "media/upper_$i.jld";
+        id_list,ts,te = JLD.load(filename,"veh_id_list","ts","te");
+        scene_list_true = replay_scenelist(f,id_list=id_list,ts=ts,te=te);
+
+        dirname = "$(name)_$(i)"
+        mkdir("lmfit/$(dirname)")
+        feat_dict = scenelist2idmfeatures(f,scene_list_true,id_list=id_list,
+                        scenario_name=name,scenario_number=i);
+end
 
 #********************Train upper test lower******************
 # We need to show a variability in the generated scenarios
@@ -113,4 +131,4 @@ coll_barchart(coll_mat_list,filename = "media/coll_barchart_upper.svg");
 #         push!(test_tdata,Dict("t"=>t))
 # end
 
-# curve_fit(test_lmfit,test_tdata,ydata,[0.5,0.5])
+# curve_fit(test_lmfit,test_tdata,ydata,[0.5,0.5])

scripts/helpers.jl

@@ -7,6 +7,7 @@ Most functions have the script within the doc string as example
 using Distributions # provides `mean`: to compute mean of particle dist over cars
 using JLD # to save models
 using AutomotiveInteraction
+using DelimitedFiles # to write to txt file that will be read in by python
 
 # We need a function to read in the jld file and extract metrics based on it
 """
@@ -92,6 +93,40 @@ modelmaker=make_cidm_models)
     return rmse_pos,rmse_vel,c_array
 end
 
+"""
+function make_lmfit_idm_models(f::FilteringEnvironment,scene)
+
+- Read lmfit params stored in txt file and create associated driver model
+
+# Examples
+cd("scripts")
+include("helpers.jl")
+f = FilteringEnvironment()
+ts = 310 # Corresponds to scenario 4
+scene = deepcopy(f.traj[ts])
+models = make_lmfit_idm_models(f,scene,scenario_name="upper",scenario_number=4)
+"""
+function make_lmfit_idm_models(f::FilteringEnvironment,scene;scenario_name="upper",
+    scenario_number=1)
+    print("Models assigned based on lmfit parameters\n")
+    models = Dict{Int64,DriverModel}()
+    for veh in scene
+        veh_id = veh.id
+        if isfile("pythonscripts/$(scenario_name)_$(scenario_number)/$(veh_id)_lmfit_params.txt")
+            lmparams = readdlm("pythonscripts/upper_4/$(veh_id)_lmfit_params.txt")
+            v_des = lmparams[1]
+            T = lmparams[2]
+            dmin = lmparams[3]
+            models[veh.id] = IntelligentDriverModel(v_des=v_des,T=T,s_min=dmin)
+        else
+            # Default to Jeremy's parameters
+            models[veh.id] = IntelligentDriverModel(v_des=17.837,s_min=5.249,
+                                T=0.918,a_max=0.758,d_cmf=3.811)
+        end
+    end
+    return models
+end
+
 """
 - Run multiple scenarios for the idm based models
 
@@ -453,29 +488,71 @@ function scenelist2idmfeatures(f,scene_list;id_list=[])
 - Returns dict with vehicle id as key and array of idm features as value
 - Array has each row a different timestep. Col 1 is v_ego, col2 is delta_v, 3 is headway
 - Col4 is the true acceleration
+- Writes the features and accleration trace to .txt files
 
 # Example
 ```julia
-f = FilteringEnvironment()
-filename = "media/upper_4.jld"
-models,id_list,ts,te = JLD.load(filename,"m","veh_id_list","ts","te")
-scene_list_true = replay_scenelist(f,id_list=id_list,ts=ts,te=te)
-feat_dict = scenelist2idmfeatures(f,scene_list_true,id_list=id_list)
+using AutomotiveInteraction
+using AutomotiveSimulator
+using JLD
+cd("scripts")
+include("helpers.jl")
+f = FilteringEnvironment();
+filename = "media/upper_4.jld";
+models,id_list,ts,te = JLD.load(filename,"m","veh_id_list","ts","te");
+scene_list_true = replay_scenelist(f,id_list=id_list,ts=ts,te=te);
+feat_dict = scenelist2idmfeatures(f,scene_list_true,id_list=id_list);
 ```
 """
-function scenelist2idmfeatures(f,scene_list;id_list=[])
+function scenelist2idmfeatures(f,scene_list;id_list=[],
+    scenario_name="upper",scenario_number=4)
     numscenes=length(scene_list)
     idmfeat_dict = Dict()
     for vehid in id_list
-        idmfeats = fill(0.,numscenes,4) # Because 3 idm features and 1 true accl
+        data = fill(0.,numscenes,4) # Because 3 idm features and 1 true accl
         acc_trace = acc(f.roadway,scene_list,vehid)
         acc_trace[1] = -1. # Just to get rid of missing. This will be thrown away later
         for (i,scene) in enumerate(scene_list)
             scene = scene_list[i]
-            idmfeats[i,1],idmfeats[i,2],idmfeats[i,3] = extract_idm_features(f,scene,vehid)
-            idmfeats[i,4] = acc_trace[i] # Caution: 1st elem will be missing
+            data[i,1],data[i,2],data[i,3] = extract_idm_features(f,scene,vehid)
+            data[i,4] = acc_trace[i] # Caution: 1st elem will be missing
         end
-        idmfeat_dict[vehid] = idmfeats
+
+        temp = removeNaNrows(data)
+        idm_feats = temp[2:end,1:3]
+        acc_trace = temp[2:end,4]
+
+        dirname = "$(scenario_name)_$(scenario_number)"
+
+        writedlm("lmfit/$(dirname)/$(vehid)_idmfeats.txt",idm_feats)
+        writedlm("lmfit/$(dirname)/$(vehid)_trueacc.txt",acc_trace)
+
+        idmfeat_dict[vehid] = temp
     end
     return idmfeat_dict
 end
+
+"""
+function removeNaNrows(a)
+- Remove rows with any NaNs in them from the input matrix `a`
+
+# Example
+```julia
+a = rand(6,3)
+a[2,2] = NaN, a[5,3]=NaN
+b = removeNaNrows(a)
+```
+"""
+function removeNaNrows(a)
+    nan_row_idx = []
+    for i in 1:size(a,1)
+        curr_row = a[i,:]
+
+        if sum(isnan.(curr_row))>0
+            push!(nan_row_idx,i)
+        end
+    end
+
+    b = a[setdiff(1:end,nan_row_idx),:]
+    return b
+end