Skip to content
Snippets Groups Projects
Commit e20a6246 authored by Patrick Godwin's avatar Patrick Godwin Committed by Kipp Cannon
Browse files

gstlal_sngltrigger: factor out code that adds single peaks to buffers to avoid... 

gstlal_sngltrigger: factor out code that adds single peaks to buffers to avoid looping through channels to select max snr channel when using this option
parent e8d2de0b
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
gstlal-burst/lib/gstlal_sngltrigger.c
+ 114
− 104
View file @ e20a6246
/*
* Copyright (C) 2012,2013,2015,2016 Chad Hanna
* 2017-2018 Duncan Meacher
*
* This program is free software; you can redistribute it and/or modify it
* under the terms of the GNU General Public License as published by the
......@@ -57,120 +58,129 @@ GstBuffer *gstlal_sngltrigger_new_buffer_from_peak(struct gstlal_peak_state *inp
GST_BUFFER_DURATION(srcbuf) = (GstClockTime) gst_util_uint64_scale_int_round(GST_SECOND, length, rate);
if (input->num_events) {
if(max_snr)
{
max_snr_index = gstlal_peak_max_over_channels(input);
/* Type casting unsigned int (guint) to int */
g_assert(max_snr_index >= 0 && max_snr_index < (int)input->channels);
}
guint channel;
for(channel = 0; channel < input->channels; channel++) {
struct GSTLALSnglTrigger *event;
SnglTriggerTable *parent;
double complex maxdata_channel = 0;
/* create a single buffer from max snr channel */
if(max_snr) {
max_snr_index = gstlal_peak_max_over_channels(input);
/* Type casting unsigned int (guint) to int */
if(max_snr && (int)channel != max_snr_index)
continue;
g_assert(max_snr_index >= 0 && max_snr_index < (int)input->channels);
switch (input->type)
{
case GSTLAL_PEAK_COMPLEX:
maxdata_channel = (double complex) input->interpvalues.as_float_complex[channel];
break;
/* create a single buffer from max snr channel */
add_buffer_from_channel(input, channel_name, pad, length, time, rate, chi2, snr_matrix_view, max_snr_index, srcbuf);
case GSTLAL_PEAK_DOUBLE_COMPLEX:
maxdata_channel = (double complex) input->interpvalues.as_double_complex[channel];
break;
/* add a buffer from each channel */
} else {
default:
g_assert(input->type == GSTLAL_PEAK_COMPLEX || input->type == GSTLAL_PEAK_DOUBLE_COMPLEX);
guint channel;
for(channel = 0; channel < input->channels; channel++) {
add_buffer_from_channel(input, channel_name, pad, length, time, rate, chi2, snr_matrix_view, channel, srcbuf);
}
}
}
if (!maxdata_channel)
continue;
/*
* allocate new event structure
*/
/*
* Populate the SNR snippet if available
* FIXME: only supported for single precision at the moment
*/
if (snr_matrix_view)
{
/* Get the column of SNR we are interested in */
gsl_vector_complex_float_view snr_vector_view = gsl_matrix_complex_float_column(&(snr_matrix_view->matrix), channel);
/* Allocate an empty time series to hold it. The event takes ownership, so no need to free it*/
event = gstlal_sngltrigger_new(snr_vector_view.vector.size);
/* Make a GSL view of the time series array data */
gsl_vector_complex_float_view snr_series_view = gsl_vector_complex_float_view_array((float *) event->snr, event->length);
/* Use BLAS to do the copy */
gsl_blas_ccopy (&(snr_vector_view.vector), &(snr_series_view.vector));
} else
event = gstlal_sngltrigger_new(0);
parent = (SnglTriggerTable *) event;
if (!event) {
/* FIXME handle error */
}
return srcbuf;
}
/*
* populate
*/
//*parent = bankarray[channel];
strcpy(parent->channel, channel_name);
parent->snr = cabs(maxdata_channel);
parent->phase = carg(maxdata_channel);
parent->channel_index = channel;
XLALINT8NSToGPS(&event->epoch, time);
{
LIGOTimeGPS end_time = event->epoch;
XLALGPSAdd(&end_time, (double) input->interpsamples[channel] / rate);
XLALGPSAddGPS(&parent->end, &end_time);
}
XLALGPSAdd(&event->epoch, (double) (input->samples[channel] - input->pad) / rate);
event->deltaT = 1. / rate;
/* populate chi squared if we have it */
parent->chisq = 0.0;
switch (input->type)
{
case GSTLAL_PEAK_COMPLEX:
if (chi2) parent->chisq = (double) *(((float *) chi2 ) + channel);
break;
case GSTLAL_PEAK_DOUBLE_COMPLEX:
if (chi2) parent->chisq = (double) *(((double *) chi2 ) + channel);
break;
default:
g_assert(input->type == GSTLAL_PEAK_COMPLEX || input->type == GSTLAL_PEAK_DOUBLE_COMPLEX);
}
void add_buffer_from_channel(struct gstlal_peak_state *input, char *channel_name, GstPad *pad, guint64 length, GstClockTime time, guint rate, void *chi2, gsl_matrix_complex_float_view *snr_matrix_view, int channel, GstBuffer *srcbuf)
{
struct GSTLALSnglTrigger *event;
SnglTriggerTable *parent;
double complex maxdata_channel = 0;
switch (input->type)
{
case GSTLAL_PEAK_COMPLEX:
maxdata_channel = (double complex) input->interpvalues.as_float_complex[channel];
break;
case GSTLAL_PEAK_DOUBLE_COMPLEX:
maxdata_channel = (double complex) input->interpvalues.as_double_complex[channel];
break;
default:
g_assert(input->type == GSTLAL_PEAK_COMPLEX || input->type == GSTLAL_PEAK_DOUBLE_COMPLEX);
}
/*
* add to buffer
*/
gst_buffer_append_memory(
srcbuf,
gst_memory_new_wrapped(
GST_MEMORY_FLAG_READONLY | GST_MEMORY_FLAG_PHYSICALLY_CONTIGUOUS,
event,
sizeof(*event) + sizeof(event->snr[0]) * event->length,
0,
sizeof(*event) + sizeof(event->snr[0]) * event->length,
event,
(GDestroyNotify) gstlal_sngltrigger_free
)
);
}
if (!maxdata_channel)
return;
/*
* allocate new event structure
*/
/*
* Populate the SNR snippet if available
* FIXME: only supported for single precision at the moment
*/
if (snr_matrix_view)
{
/* Get the column of SNR we are interested in */
gsl_vector_complex_float_view snr_vector_view = gsl_matrix_complex_float_column(&(snr_matrix_view->matrix), channel);
/* Allocate an empty time series to hold it. The event takes ownership, so no need to free it*/
event = gstlal_sngltrigger_new(snr_vector_view.vector.size);
/* Make a GSL view of the time series array data */
gsl_vector_complex_float_view snr_series_view = gsl_vector_complex_float_view_array((float *) event->snr, event->length);
/* Use BLAS to do the copy */
gsl_blas_ccopy (&(snr_vector_view.vector), &(snr_series_view.vector));
} else
event = gstlal_sngltrigger_new(0);
parent = (SnglTriggerTable *) event;
if (!event) {
/* FIXME handle error */
}
return srcbuf;
/*
* populate
*/
//*parent = bankarray[channel];
strcpy(parent->channel, channel_name);
parent->snr = cabs(maxdata_channel);
parent->phase = carg(maxdata_channel);
parent->channel_index = channel;
XLALINT8NSToGPS(&event->epoch, time);
{
LIGOTimeGPS end_time = event->epoch;
XLALGPSAdd(&end_time, (double) input->interpsamples[channel] / rate);
XLALGPSAddGPS(&parent->end, &end_time);
}
XLALGPSAdd(&event->epoch, (double) (input->samples[channel] - input->pad) / rate);
event->deltaT = 1. / rate;
/* populate chi squared if we have it */
parent->chisq = 0.0;
switch (input->type)
{
case GSTLAL_PEAK_COMPLEX:
if (chi2) parent->chisq = (double) *(((float *) chi2 ) + channel);
break;
case GSTLAL_PEAK_DOUBLE_COMPLEX:
if (chi2) parent->chisq = (double) *(((double *) chi2 ) + channel);
break;
default:
g_assert(input->type == GSTLAL_PEAK_COMPLEX || input->type == GSTLAL_PEAK_DOUBLE_COMPLEX);
}
/*
* add to buffer
*/
gst_buffer_append_memory(
srcbuf,
gst_memory_new_wrapped(
GST_MEMORY_FLAG_READONLY | GST_MEMORY_FLAG_PHYSICALLY_CONTIGUOUS,
event,
sizeof(*event) + sizeof(event->snr[0]) * event->length,
0,
sizeof(*event) + sizeof(event->snr[0]) * event->length,
event,
(GDestroyNotify) gstlal_sngltrigger_free
)
);
}
gstlal-burst/lib/gstlal_sngltrigger.h
+ 1
− 0
View file @ e20a6246
......@@ -65,6 +65,7 @@ SnglTriggerTable;
*/
GstBuffer *gstlal_sngltrigger_new_buffer_from_peak(struct gstlal_peak_state *input, char *channel_name, GstPad *pad, guint64 offset, guint64 length, GstClockTime time, guint rate, void *chi2, gsl_matrix_complex_float_view *snr_matrix_view, GstClockTimeDiff, gboolean max_snr);
void add_buffer_from_channel(struct gstlal_peak_state *input, char *channel_name, GstPad *pad, guint64 length, GstClockTime time, guint rate, void *chi2, gsl_matrix_complex_float_view *snr_matrix_view, int channel, GstBuffer *srcbuf);
G_END_DECLS
#endif /* __GSTLAL_SNGLTRIGGER_H__ */
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment