Some work on the mlat stuff, and restoring sanity to the gr-uhd tags interface
This commit is contained in:
@@ -80,14 +80,15 @@ static double correlate_preamble(const float *in, int samples_per_chip) {
|
||||
|
||||
|
||||
static double pmt_to_timestamp(pmt::pmt_t tstamp, uint64_t abs_sample_cnt, double secs_per_sample) {
|
||||
uint64_t ts_sample, ticks;
|
||||
uint64_t ts_sample;
|
||||
double last_stamp;
|
||||
|
||||
if(pmt::pmt_symbol_to_string(gr_tags::get_key(tstamp)) != "timestamp") return 0;
|
||||
|
||||
ticks = pmt_to_uint64(gr_tags::get_value(tstamp));
|
||||
last_stamp = pmt_to_double(gr_tags::get_value(tstamp));
|
||||
ts_sample = gr_tags::get_nitems(tstamp);
|
||||
//std::cout << "HEY WE GOT A STAMP AT " << ticks << " TICKS AT SAMPLE " << ts_sample << " ABS SAMPLE CNT IS " << abs_sample_cnt << std::endl;
|
||||
return double(ticks+abs_sample_cnt) * secs_per_sample;
|
||||
return double(abs_sample_cnt * secs_per_sample) + last_stamp;
|
||||
}
|
||||
|
||||
int air_modes_preamble::general_work(int noutput_items,
|
||||
@@ -114,7 +115,7 @@ int air_modes_preamble::general_work(int noutput_items,
|
||||
if(tstamp_tags.size() > 0) {
|
||||
d_timestamp = tstamp_tags.back();
|
||||
}
|
||||
|
||||
|
||||
for(int i=0; i < ninputs; i++) {
|
||||
float pulse_threshold = inavg[i] * d_threshold;
|
||||
if(in[i] > pulse_threshold) { //hey we got a candidate
|
||||
@@ -125,15 +126,18 @@ int air_modes_preamble::general_work(int noutput_items,
|
||||
if( in[i+pulse_offsets[3]] < pulse_threshold ) continue;
|
||||
|
||||
//get a more accurate bit center by finding the correlation peak across all four preamble bits
|
||||
bool late = false;
|
||||
bool late, early;
|
||||
do {
|
||||
double now_corr = correlate_preamble(in+i, d_samples_per_chip);
|
||||
double late_corr = correlate_preamble(in+i+1, d_samples_per_chip);
|
||||
double early_corr = correlate_preamble(in+i-1, d_samples_per_chip);
|
||||
late = (late_corr > now_corr);
|
||||
//early = (early_corr > now_corr);
|
||||
if(late) i++;
|
||||
} while(late);
|
||||
//if(early && i>0) { std::cout << "EARLY " << i << std::endl; i--; }
|
||||
} while(late);// xor early);
|
||||
|
||||
//now check to see that the rest of the chips in the preamble
|
||||
//now check to see that the non-peak symbols in the preamble
|
||||
//are below the peaks by threshold dB
|
||||
float avgpeak = ( in[i+pulse_offsets[0]]
|
||||
+ in[i+pulse_offsets[1]]
|
||||
|
||||
@@ -79,16 +79,36 @@ for msg in position_reports:
|
||||
offset_list.append({"aircraft": data["shortdata"], "times": rel_times})
|
||||
|
||||
#this is a list of unique aircraft, heard by all stations, which transmitted position packets
|
||||
#we do drift calcs separately for each aircraft in the set because mixing them seems to screw things up
|
||||
#i haven't really sat down and figured out why that is yet
|
||||
unique_aircraft = list(set([x["aircraft"] for x in offset_list]))
|
||||
#todo: the below can be done cleaner with nested list comprehensions
|
||||
for ac in unique_aircraft:
|
||||
for i in range(1,len(stations)):
|
||||
#pull out a list of unique aircraft from the offset list
|
||||
rel_times_for_one_ac = [report["times"][i]-report["times"][0] for report in offset_list if report["aircraft"] == ac]
|
||||
abs_times_for_one_ac = [report["times"][0] for report in offset_list if report["aircraft"] == ac]
|
||||
|
||||
#find drift error
|
||||
drift_error = [(y-x)/(b-a) for x,y,a,b in zip(rel_times_for_one_ac, rel_times_for_one_ac[1:], abs_times_for_one_ac, abs_times_for_one_ac[1:])]
|
||||
drift_error_limited = [x for x in drift_error if abs(x) < 1e-5]
|
||||
print "drift from %d relative to station 0 for ac %x: %.3fppm" % (i, ac & 0xFFFFFF, numpy.mean(drift_error_limited) * 1e6)
|
||||
|
||||
#get a list of reported times gathered by the unique aircraft that transmitted them
|
||||
#abs_unique_times = [report["times"] for ac in unique_aircraft for report in offset_list if report["aircraft"] == ac]
|
||||
#print abs_unique_times
|
||||
#todo: the below can be done cleaner with nested list comprehensions
|
||||
clock_rate_corrections = [0]
|
||||
for i in range(1,len(stations)):
|
||||
drift_error_limited = []
|
||||
for ac in unique_aircraft:
|
||||
times = [report["times"] for report in offset_list if report["aircraft"] == ac]
|
||||
|
||||
s0_times = [report[0] for report in times]
|
||||
rel_times = [report[i]-report[0] for report in times]
|
||||
|
||||
#find drift error rate
|
||||
drift_error = [(y-x)/(b-a) for x,y,a,b in zip(rel_times, rel_times[1:], s0_times[0:], s0_times[1:])]
|
||||
drift_error_limited.append([x for x in drift_error if abs(x) < 1e-5])
|
||||
|
||||
#flatten the list of lists (tacky, there's a better way)
|
||||
drift_error_limited = [x for sublist in drift_error_limited for x in sublist]
|
||||
clock_rate_corrections.append(0-numpy.mean(drift_error_limited))
|
||||
|
||||
for i in range(len(clock_rate_corrections)):
|
||||
print "drift from %d relative to station 0: %.3fppm" % (i, clock_rate_corrections[i] * 1e6)
|
||||
|
||||
#ok now we have clock rate corrections for each station, and we can multilaterate like hell
|
||||
#we'll need a way to get a clock offset at a particular time, in case the dataset has a discontinuity
|
||||
#actually we don't, unless we still have the old "timestamps-are-wrong-in-b100" issue
|
||||
#remember to subtract out rel_time before correcting for clock rate error
|
||||
#start with all_heard
|
||||
|
||||
Reference in New Issue
Block a user