cycling/ros_variance.m

   1 function [ros,ros_mean] = ros_variance(n,ps)
   2 %estimates variance in ros between two detections
   3 %n - number of random points to use in simulation
   4 %ps - struct with path information, ps = cluster_paths(w,1)
   5
   6 real_data = input_num('Use real data? 0 = no 1 = yes',1,1);
   7
   8 if real_data == 0
   9     %coordinates of detections
  10     d1 = [0 0 0];
  11     d2 = [3 4 12];
  12
  13     %reported
  14     dist = sqrt((d2(1)-d1(1))^2+(d2(2)-d1(2))^2);
  15     time = d2(3) - d1(3);
  16     ros = dist/time;
  17     fprintf('ROS = %f\n',ros)
  18
  19     %statistics
  20     sig1 = 1;
  21     sig2 = 1;
  22     time_diff = 0;
  23
  24 else
  25     %real data
  26     %load p_struct.mat
  27     %ps = p;
  28     %clear p
  29     num_paths = length(ps.paths);
  30     path_num = max(2,round(num_paths*rand));
  31     path = ps.paths(path_num);
  32     path_length = length(path);
  33     path_end = max(2,round(path_length*rand));
  34     p1 = path.p(path_end -1);
  35     p2 = path.p(path_end);
  36     fprintf('Generating stats for path %d, points %d and %d \n',path_num,path_end-1,path_end);
  37     lon1 = ps.points(p1,2);
  38     lat1 = ps.points(p1,1);
  39     time1 = ps.points(p1,3);
  40     lon2 = ps.points(p2,2);
  41     lat2 = ps.points(p2,1);
  42     time2 = ps.points(p2,3);
  43
  44     E = wgs84Ellipsoid;
  45     dist = distance(lat1,lon1,lat2,lon2,E);
  46     time = time2-time1;
  47     ros = dist/(time*24*3600);
  48     fprintf('ROS = %f\n',ros)
  49
  50     %coordinates of detections
  51     d1 = [0 0 0];
  52     d2 = [0 dist time];
  53
  54     %statistics
  55     sig1 = 1000;
  56     sig2 = 1000;
  57     time_diff = 0.25;
  58 end %if real_data == 0
  59
  60
  61
  62 %n = 200;
  63 r1 = sig1*randn(n,3)+d1;
  64 r1(:,3) = d1(3)-time_diff*rand(n,1);
  65 r2 = sig2*randn(n,3)+d2;
  66 r2(:,3) = d2(3)-time_diff*rand(n,1);
  67
  68 dr = r2-r1;
  69 dist_rand = sqrt(dr(:,1).^2+dr(:,2).^2);
  70 time_rand = dr(:,3);
  71 if real_data == 0
  72     ros_rand = dist_rand./time_rand;
  73 else
  74     ros_rand = dist_rand./(24*3600*time_rand);
  75 end
  76 if n < 1000
  77     figure,scatter3(r1(:,1),r1(:,2),r1(:,3),'*b')
  78     xlabel('x'),ylabel('y'),zlabel('time')
  79     title('Random paths between two detections')
  80     hold on
  81     scatter3(r2(:,1),r2(:,2),r2(:,3),'*r')
  82     for i = 1:n
  83         plot3([r1(i,1) r2(i,1)],[r1(i,2) r2(i,2)],[r1(i,3) r2(i,3)],':')
  84     end
  85     hold off
  86     %make joint pdf
  87     figure,scatter3(dist_rand,1./time_rand,ros_rand);
  88     xlabel('Distance'),ylabel('1/Time'),zlabel('ROS')
  89 end
  90
  91 figure,histogram(ros_rand)
  92 if real_data == 0
  93     title('Histogram of ROS')
  94 else
  95     t_str = sprintf('Histogram of ROS \n Path: %d Points: %d and %d',path_num,path_end-1,path_end);
  96     title(t_str)
  97 end
  98 xlabel('ROS'),ylabel('Number')
  99 ros_mean = mean(ros_rand);
 100 fprintf('Mean random ROS: %f\n',mean(ros_rand))
 101
 102
 103 end