diff --git a/Appl/CEP/CS1/CS1_BGLProc/src/FIR.S b/Appl/CEP/CS1/CS1_BGLProc/src/FIR.S index 48f53e02ac7daa9c5a757881246dabb86267b311..cadf5c0838519b0aa95ac4a180d696e96e52d736 100644 --- a/Appl/CEP/CS1/CS1_BGLProc/src/FIR.S +++ b/Appl/CEP/CS1/CS1_BGLProc/src/FIR.S @@ -32,6 +32,7 @@ _FIR_constants_used: .long NR_STATIONS .long NR_SAMPLES_PER_INTEGRATION .long NR_SUBBAND_CHANNELS + .long NR_TAPS .long NR_POLARIZATIONS diff --git a/Appl/CEP/CS1/CS1_BGLProc/src/FIR.h b/Appl/CEP/CS1/CS1_BGLProc/src/FIR.h index 12640c2dba05fe55903b23c462db8abd30beb687..08db234efff2e806bda07432bcc3afecf4298d7e 100644 --- a/Appl/CEP/CS1/CS1_BGLProc/src/FIR.h +++ b/Appl/CEP/CS1/CS1_BGLProc/src/FIR.h @@ -69,6 +69,7 @@ extern "C" { unsigned nr_stations; unsigned nr_samples_per_integration; unsigned nr_subband_channels; + unsigned nr_taps; unsigned nr_polarizations; } _FIR_constants_used; diff --git a/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.cc b/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.cc index e2166636455a94c2450c18308b6e5d8db301d7f9..ca669e6c670b8abf2dccb7bc605753a6c2a80d27 100644 --- a/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.cc +++ b/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.cc @@ -41,10 +41,8 @@ #include <rts.h> #endif -namespace LOFAR -{ -namespace CS1 -{ +namespace LOFAR { +namespace CS1 { #if !defined HAVE_MASS @@ -57,7 +55,11 @@ inline static dcomplex cosisin(double x) FIR WH_BGL_Processing::itsFIRs[NR_STATIONS][NR_POLARIZATIONS][NR_SUBBAND_CHANNELS] CACHE_ALIGNED; fcomplex WH_BGL_Processing::samples[NR_SUBBAND_CHANNELS][NR_STATIONS][NR_SAMPLES_PER_INTEGRATION][NR_POLARIZATIONS] CACHE_ALIGNED; +#if defined SPARSE_FLAGS +SparseSet WH_BGL_Processing::flags[NR_STATIONS]; +#else bitset<NR_SAMPLES_PER_INTEGRATION> WH_BGL_Processing::flags[NR_STATIONS] CACHE_ALIGNED; +#endif unsigned WH_BGL_Processing::itsNrValidSamples[NR_BASELINES] CACHE_ALIGNED; float WH_BGL_Processing::correlationWeights[NR_SAMPLES_PER_INTEGRATION + 1] CACHE_ALIGNED; float WH_BGL_Processing::thresholds[NR_BASELINES][NR_SUBBAND_CHANNELS]; @@ -1145,6 +1147,18 @@ WH_BGL_Processing::WH_BGL_Processing(const string& name, double baseFrequency, c #if !defined C_IMPLEMENTATION ASSERT(NR_SAMPLES_PER_INTEGRATION % 16 == 0); + + ASSERT(_FIR_constants_used.input_type == INPUT_TYPE); + ASSERT(_FIR_constants_used.nr_stations == NR_STATIONS); + ASSERT(_FIR_constants_used.nr_samples_per_integration == NR_SAMPLES_PER_INTEGRATION); + ASSERT(_FIR_constants_used.nr_subband_channels == NR_SUBBAND_CHANNELS); + ASSERT(_FIR_constants_used.nr_taps == NR_TAPS); + ASSERT(_FIR_constants_used.nr_polarizations == NR_POLARIZATIONS); + + ASSERT(_correlator_constants_used.nr_stations == NR_STATIONS); + ASSERT(_correlator_constants_used.nr_samples_per_integration == NR_SAMPLES_PER_INTEGRATION); + ASSERT(_correlator_constants_used.nr_subband_channels == NR_SUBBAND_CHANNELS); + ASSERT(_correlator_constants_used.nr_polarizations == NR_POLARIZATIONS); #endif getDataManager().addInDataHolder(SUBBAND_CHANNEL, new DH_Subband("input", ps)); @@ -1194,17 +1208,6 @@ void WH_BGL_Processing::preprocess() } #if defined HAVE_BGL && !defined C_IMPLEMENTATION - ASSERT(_FIR_constants_used.input_type == INPUT_TYPE); - ASSERT(_FIR_constants_used.nr_stations == NR_STATIONS); - ASSERT(_FIR_constants_used.nr_samples_per_integration == NR_SAMPLES_PER_INTEGRATION); - ASSERT(_FIR_constants_used.nr_subband_channels == NR_SUBBAND_CHANNELS); - ASSERT(_FIR_constants_used.nr_polarizations == NR_POLARIZATIONS); - - ASSERT(_correlator_constants_used.nr_stations == NR_STATIONS); - ASSERT(_correlator_constants_used.nr_samples_per_integration == NR_SAMPLES_PER_INTEGRATION); - ASSERT(_correlator_constants_used.nr_subband_channels == NR_SUBBAND_CHANNELS); - ASSERT(_correlator_constants_used.nr_polarizations == NR_POLARIZATIONS); - mutex = rts_allocate_mutex(); #endif } @@ -1212,17 +1215,32 @@ void WH_BGL_Processing::preprocess() void WH_BGL_Processing::computeFlags() { + computeFlagsTimer.start(); + #if NR_SUBBAND_CHANNELS == 1 DH_Subband::AllFlagsType *input = get_DH_Subband()->getFlags(); bitset<NR_INPUT_SAMPLES> (&flags)[NR_STATIONS] = *input; +#else +#if defined SPARSE_FLAGS + DH_Subband::AllFlagsType *input = get_DH_Subband()->getFlags(); + + for (int stat = 0; stat < NR_STATIONS; stat ++) { + flags[stat].reset(); + const std::vector<SparseSet::range> &ranges = (*input)[stat].getRanges(); + + for (std::vector<SparseSet::range>::const_iterator it = ranges.begin(); it != ranges.end(); it ++) { + int begin = std::max(0, (int) it->begin / NR_SUBBAND_CHANNELS - NR_TAPS + 1); + int end = std::min(NR_SAMPLES_PER_INTEGRATION - 1, (int) (it->end - 1) / NR_SUBBAND_CHANNELS); + + flags[stat].include(begin, end); + } + } #else typedef bitset<NR_SUBBAND_CHANNELS> inputType[NR_STATIONS][NR_TAPS - 1 + NR_SAMPLES_PER_INTEGRATION]; inputType *input = (inputType *) get_DH_Subband()->getFlags(); - computeFlagsTimer.start(); - memset(flags, 0, sizeof flags); for (int stat = 0; stat < NR_STATIONS; stat ++) { @@ -1244,6 +1262,7 @@ void WH_BGL_Processing::computeFlags() } } } +#endif #endif for (int stat2 = 0; stat2 < NR_STATIONS; stat2 ++) { @@ -1330,7 +1349,11 @@ void WH_BGL_Processing::doPPF() FFTtimer.start(); for (int time = 0; time < NR_SAMPLES_PER_INTEGRATION; time ++) { for (int pol = 0; pol < NR_POLARIZATIONS; pol ++) { +#if defined SPARSE_FLAGS + if (flags[stat].test(time)) { +#else if (flags[stat][time]) { +#endif for (int chan = 0; chan < NR_SUBBAND_CHANNELS; chan ++) { samples[chan][stat][time][pol] = makefcomplex(0, 0); } @@ -1380,6 +1403,67 @@ void WH_BGL_Processing::doPPF() computePhaseShifts(phaseShifts, (*delays)[stat]); #endif +#if defined SPARSE_FLAGS + const std::vector<SparseSet::range> &ranges = flags[stat].getRanges(); + std::vector<SparseSet::range>::const_iterator it = ranges.begin(); + + for (int time = 0; time < NR_SAMPLES_PER_INTEGRATION; it ++) { + int firstBad, firstGood; + + if (it != ranges.end()) { + firstBad = it->begin; + firstGood = it->end; + } else { + firstBad = firstGood = NR_SAMPLES_PER_INTEGRATION; + } + + for (; time < firstBad; time ++) { + FFTtimer.start(); + _prefetch(fftInData[time], + sizeof(fcomplex[NR_POLARIZATIONS][NR_SUBBAND_CHANNELS]) / CACHE_LINE_SIZE, + CACHE_LINE_SIZE); + + for (int pol = 0; pol < NR_POLARIZATIONS; pol ++) { + fftw_one(itsFFTWPlan, + (fftw_complex *) fftInData[time][pol], + (fftw_complex *) fftOutData[time & 1][pol]); + } + FFTtimer.stop(); + + if (time & 1) { +#if defined DELAY_COMPENSATION + _phase_shift_and_transpose(&samples[0][stat][time - 1][0], + &fftOutData[0][0][0], + &phaseShifts[time - 1]); +#else + _transpose_4x8(&samples[0][stat][time - 1][0], + &fftOutData[0][0][0], + NR_SUBBAND_CHANNELS, + sizeof(fcomplex) * NR_SUBBAND_CHANNELS, + sizeof(fcomplex) * NR_POLARIZATIONS * NR_SAMPLES_PER_INTEGRATION * NR_STATIONS); +#endif + } + } + + for (; time < firstGood; time ++) { + _memzero(fftOutData[time & 1], sizeof(fftOutData[time & 1])); + + if (time & 1) { +#if defined DELAY_COMPENSATION + _phase_shift_and_transpose(&samples[0][stat][time - 1][0], + &fftOutData[0][0][0], + &phaseShifts[time - 1]); +#else + _transpose_4x8(&samples[0][stat][time - 1][0], + &fftOutData[0][0][0], + NR_SUBBAND_CHANNELS, + sizeof(fcomplex) * NR_SUBBAND_CHANNELS, + sizeof(fcomplex) * NR_POLARIZATIONS * NR_SAMPLES_PER_INTEGRATION * NR_STATIONS); +#endif + } + } + } +#else // !SPARSE_FLAGS for (int time = 0; time < NR_SAMPLES_PER_INTEGRATION; time += 2) { FFTtimer.start(); for (int t = 0; t < 2; t ++) { @@ -1411,7 +1495,8 @@ void WH_BGL_Processing::doPPF() sizeof(fcomplex) * NR_POLARIZATIONS * NR_SAMPLES_PER_INTEGRATION * NR_STATIONS); #endif } -#endif +#endif // SPARSE_FLAGS +#endif // C_IMPLEMENTATION } #if defined HAVE_BGL && !defined C_IMPLEMENTATION @@ -1632,6 +1717,11 @@ void WH_BGL_Processing::doCorrelate() void WH_BGL_Processing::process() { totalTimer.start(); + +#if defined SPARSE_FLAGS + get_DH_Subband()->getExtraData(); +#endif + computeFlags(); #if NR_SUBBAND_CHANNELS > 1 @@ -1655,6 +1745,5 @@ void WH_BGL_Processing::dump() const { } - } // namespace CS1 } // namespace LOFAR diff --git a/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.h b/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.h index 55821a7d096fc1176eca6b00512e11cb23c499ae..7759e1d30055d8edbd6262843282d86c40340e8e 100644 --- a/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.h +++ b/Appl/CEP/CS1/CS1_BGLProc/src/WH_BGL_Processing.h @@ -127,7 +127,11 @@ class WH_BGL_Processing: public WorkHolder { static FIR itsFIRs[NR_STATIONS][NR_POLARIZATIONS][NR_SUBBAND_CHANNELS] CACHE_ALIGNED; static fcomplex samples[NR_SUBBAND_CHANNELS][NR_STATIONS][NR_SAMPLES_PER_INTEGRATION][NR_POLARIZATIONS] CACHE_ALIGNED; +#if defined SPARSE_FLAGS + static SparseSet WH_BGL_Processing::flags[NR_STATIONS]; +#else static bitset<NR_SAMPLES_PER_INTEGRATION> flags[NR_STATIONS] CACHE_ALIGNED; +#endif static unsigned itsNrValidSamples[NR_BASELINES] CACHE_ALIGNED; static float correlationWeights[NR_SAMPLES_PER_INTEGRATION + 1] CACHE_ALIGNED; static float thresholds[NR_BASELINES][NR_SUBBAND_CHANNELS]; diff --git a/Appl/CEP/CS1/CS1_BGLProc/test/tWH_BGL_Processing.cc b/Appl/CEP/CS1/CS1_BGLProc/test/tWH_BGL_Processing.cc index 4337cade99adf53655fe657d1e9f83aba7118570..9878986f8120544dc87e9625475db3bf7e32f369 100644 --- a/Appl/CEP/CS1/CS1_BGLProc/test/tWH_BGL_Processing.cc +++ b/Appl/CEP/CS1/CS1_BGLProc/test/tWH_BGL_Processing.cc @@ -105,11 +105,23 @@ void setSubbandTestPattern(WH_BGL_Processing &wh, double signalFrequency) #endif } +#if defined SPARSE_FLAGS + for (int stat = 0; stat < NR_STATIONS; stat ++) { + (*flags)[stat].reset(); + } +#else memset(flags, 0, sizeof(DH_Subband::AllFlagsType)); +#endif -#if 0 && NR_INPUT_SAMPLES >= 17000 +#if 1 && NR_INPUT_SAMPLES >= 17000 +#if defined SPARSE_FLAGS + //(*flags)[4].include(14000); + //(*flags)[5].include(17000); + dh->fillExtraData(); +#else (*flags)[4][14000] = true; (*flags)[5][17000] = true; +#endif #endif (std::cerr << "done.\n").flush(); diff --git a/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/DH_Subband.h b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/DH_Subband.h index e85a873dff8207424fe2d82a107c154d7492ae0d..7ecd1c0e6fdac2bfd64ae1320183589cc47d7fe0 100644 --- a/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/DH_Subband.h +++ b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/DH_Subband.h @@ -23,14 +23,19 @@ #ifndef LOFAR_CS1_INTERFACE_DH_SUBBAND_H #define LOFAR_CS1_INTERFACE_DH_SUBBAND_H +#define SPARSE_FLAGS + #include <CS1_Interface/CS1_Config.h> #include <CS1_Interface/RectMatrix.h> +#if defined SPARSE_FLAGS +#include <CS1_Interface/SparseSet.h> +#else #include <CS1_Interface/bitset.h> +#endif #include <Transport/DataHolder.h> #include <Common/lofar_complex.h> #include <APS/ParameterSet.h> - namespace LOFAR { namespace CS1 { @@ -83,10 +88,12 @@ class DH_Subband: public DataHolder return itsDelays[station]; } +#if !defined SPARSE_FLAGS const size_t nrFlags() const { return itsNrStations * ((itsNrInputSamples + 31) & ~31); } +#endif const size_t nrDelays() const { @@ -98,7 +105,11 @@ class DH_Subband: public DataHolder typedef SampleType AllSamplesType[NR_STATIONS][NR_INPUT_SAMPLES][NR_POLARIZATIONS]; // Flags +#if defined SPARSE_FLAGS + typedef SparseSet AllFlagsType[NR_STATIONS]; +#else typedef bitset<NR_INPUT_SAMPLES> AllFlagsType[NR_STATIONS]; +#endif // Fine-grained delays typedef DelayIntervalType AllDelaysType[NR_STATIONS]; @@ -136,17 +147,24 @@ class DH_Subband: public DataHolder void swapBytes(); +#if defined SPARSE_FLAGS + void getExtraData(), fillExtraData(); +#endif + private: /// Forbid assignment. DH_Subband &operator = (const DH_Subband &); - unsigned itsNrStations; - unsigned itsNrInputSamples; + unsigned itsNrStations, itsNrInputSamples; SampleType *itsSamples; // RectMatrix cannot be used for bitsets, thus not for flags RectMatrix<SampleType> *itsSamplesMatrix; +#if defined SPARSE_FLAGS + SparseSet *itsFlags; +#else uint32 *itsFlags; +#endif DelayIntervalType *itsDelays; void fillDataPointers(); diff --git a/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/Makefile.am b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/Makefile.am index 571db025d59605094bc06bb2d240ac4df1785843..9d27ce16e00f56632492c63dcc2c3abdca9d9ea3 100644 --- a/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/Makefile.am +++ b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/Makefile.am @@ -10,6 +10,7 @@ pkginclude_HEADERS = \ RectMatrix.h \ RSPTimeStamp.h \ RectMatrix.tcc \ + SparseSet.h \ Stub_BGL.h \ Stub_Delay.h \ Stub_BGL_RFI_Mitigation.h \ diff --git a/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/SparseSet.h b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/SparseSet.h new file mode 100644 index 0000000000000000000000000000000000000000..d5f4e2bcfb4acfd87e96c1d4405b078f9285eebb --- /dev/null +++ b/Appl/CEP/CS1/CS1_Interface/include/CS1_Interface/SparseSet.h @@ -0,0 +1,84 @@ +//# SparseSet.h: portable <bitset> adaptation +//# +//# Copyright (C) 2006 +//# ASTRON (Netherlands Foundation for Research in Astronomy) +//# P.O.Box 2, 7990 AA Dwingeloo, The Netherlands, seg@astron.nl +//# +//# 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 Free Software Foundation; either version 2 of the License, or +//# (at your option) any later version. +//# +//# This program is distributed in the hope that it will be useful, +//# but WITHOUT ANY WARRANTY; without even the implied warranty of +//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +//# GNU General Public License for more details. +//# +//# You should have received a copy of the GNU General Public License +//# along with this program; if not, write to the Free Software +//# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +//# +//# $Id$ + + +#ifndef LOFAR_APPL_CEP_CS1_CS1_INTERFACE_BITSET_H +#define LOFAR_APPL_CEP_CS1_CS1_INTERFACE_BITSET_H + +#include <Blob/BlobOStream.h> +#include <Blob/BlobIStream.h> + +#include <vector> + + +namespace LOFAR +{ + +class SparseSet { + public: + struct range { + unsigned begin, end; + }; + + SparseSet &include(unsigned index); + SparseSet &include(unsigned first, unsigned last); + + void reset(); + + unsigned count() const; + bool any(unsigned first, unsigned last) const; + bool test(unsigned index) const; + + SparseSet operator | (const SparseSet &) const; + SparseSet &operator -= (size_t count); + //SparseSet subset(unsigned first, unsigned last) const; /* not tested */ + + const std::vector<struct range> &getRanges() const; + + void write(BlobOStream &) const; + void read(BlobIStream &); + + private: + std::vector<struct range> ranges; +}; + + +inline SparseSet &SparseSet::include(unsigned index) +{ + return include(index, index); +} + + +inline void SparseSet::reset() +{ + ranges.resize(0); +} + + +inline const std::vector<struct SparseSet::range> &SparseSet::getRanges() const +{ + return ranges; +} + +} + +#endif diff --git a/Appl/CEP/CS1/CS1_Interface/src/DH_Subband.cc b/Appl/CEP/CS1/CS1_Interface/src/DH_Subband.cc index 09c076c6d3e55dd16bd36a9da2920c63863bcb51..a502d7d17acf8d6fb9fbeb18c7e94b9911e945ca 100644 --- a/Appl/CEP/CS1/CS1_Interface/src/DH_Subband.cc +++ b/Appl/CEP/CS1/CS1_Interface/src/DH_Subband.cc @@ -24,75 +24,117 @@ #include <Common/DataConvert.h> #include <Common/Timer.h> -namespace LOFAR +#if defined SPARSE_FLAGS +#include <Blob/BlobOStream.h> +#include <Blob/BlobIStream.h> +#endif + + +namespace LOFAR { +namespace CS1 { + +DH_Subband::DH_Subband(const string &name, const ACC::APS::ParameterSet &pSet) + : DataHolder(name, "DH_Subband"), + itsNrStations(pSet.getUint32("Observation.NStations")), + itsNrInputSamples(pSet.getUint32("Observation.NSubbandSamples") + (pSet.getUint32("BGLProc.NPPFTaps") - 1) * pSet.getUint32("Observation.NChannels")), + itsSamples(0), + itsSamplesMatrix(0), + itsFlags(0), + itsDelays(0) { - namespace CS1 - { - - DH_Subband::DH_Subband(const string &name, const ACC::APS::ParameterSet &pSet) - : DataHolder(name, "DH_Subband"), - itsNrStations(pSet.getUint32("Observation.NStations")), - itsNrInputSamples(pSet.getUint32("Observation.NSubbandSamples") + (pSet.getUint32("BGLProc.NPPFTaps") - 1) * pSet.getUint32("Observation.NChannels")), - itsSamples(0), - itsSamplesMatrix(0), - itsFlags(0), - itsDelays(0) - { - } - - DH_Subband::DH_Subband(const DH_Subband &that) - : DataHolder(that), - itsNrStations(that.itsNrStations), - itsNrInputSamples(that.itsNrInputSamples), - itsSamples(that.itsSamples), - itsSamplesMatrix(that.itsSamplesMatrix), - itsFlags(that.itsFlags), - itsDelays(that.itsDelays) - { - } - - DH_Subband::~DH_Subband() - { - delete itsSamplesMatrix; - } - - DataHolder *DH_Subband::clone() const - { - return new DH_Subband(*this); - } - - void DH_Subband::init() - { - addField("Samples", BlobField<uint8>(1, nrSamples() * sizeof(SampleType)), 32); - addField("Flags", BlobField<uint32>(1, nrFlags() / sizeof(uint32))); - addField("Delays", BlobField<float>(1, nrDelays() * sizeof(DelayIntervalType) / sizeof(float))); - - createDataBlock(); - - vector<DimDef> vdd; - vdd.push_back(DimDef("Station", itsNrStations)); - vdd.push_back(DimDef("Time", itsNrInputSamples)); - vdd.push_back(DimDef("Polarisation", NR_POLARIZATIONS)); - - itsSamplesMatrix = new RectMatrix<SampleType> (vdd); - itsSamplesMatrix->setBuffer(itsSamples, nrSamples()); - - memset(itsFlags, 0, sizeof *itsFlags); - } - - void DH_Subband::fillDataPointers() - { - itsSamples = (SampleType *) getData<uint8> ("Samples"); - itsFlags = (uint32 *) getData<uint32>("Flags"); - itsDelays = (DelayIntervalType *) getData<float> ("Delays"); - } - - void DH_Subband::swapBytes() - { - // only convert Samples; CEPframe converts Flags and Delays - dataConvert(LittleEndian, itsSamples, nrSamples()); - } - - } // namespace CS1 +#if defined SPARSE_FLAGS + setExtraBlob("Flags", 0); +#endif +} + +DH_Subband::DH_Subband(const DH_Subband &that) + : DataHolder(that), + itsNrStations(that.itsNrStations), + itsNrInputSamples(that.itsNrInputSamples), + itsSamples(that.itsSamples), + itsSamplesMatrix(that.itsSamplesMatrix), + itsFlags(that.itsFlags), + itsDelays(that.itsDelays) +{ +#if defined SPARSE_FLAGS + setExtraBlob("Flags", 0); +#endif +} + +DH_Subband::~DH_Subband() +{ +#if defined SPARSE_FLAGS + delete [] itsFlags; +#endif + delete itsSamplesMatrix; +} + +DataHolder *DH_Subband::clone() const +{ + return new DH_Subband(*this); +} + +void DH_Subband::init() +{ + addField("Samples", BlobField<uint8>(1, nrSamples() * sizeof(SampleType)), 32); + addField("Delays", BlobField<float>(1, nrDelays() * sizeof(DelayIntervalType) / sizeof(float))); +#if !defined SPARSE_FLAGS + addField("Flags", BlobField<uint32>(1, nrFlags() / sizeof(uint32))); +#endif + + createDataBlock(); + + vector<DimDef> vdd; + vdd.push_back(DimDef("Station", itsNrStations)); + vdd.push_back(DimDef("Time", itsNrInputSamples)); + vdd.push_back(DimDef("Polarisation", NR_POLARIZATIONS)); + + itsSamplesMatrix = new RectMatrix<SampleType> (vdd); + itsSamplesMatrix->setBuffer(itsSamples, nrSamples()); + +#if defined SPARSE_FLAGS + itsFlags = new SparseSet[itsNrStations]; + std::cerr << "new SparseSet[" << itsNrStations << "] = " << itsFlags << '\n'; +#endif +} + +void DH_Subband::fillDataPointers() +{ + itsSamples = (SampleType *) getData<uint8> ("Samples"); + itsDelays = (DelayIntervalType *) getData<float> ("Delays"); +#if !defined SPARSE_FLAGS + itsFlags = (uint32 *) getData<uint32>("Flags"); +#endif +} + + +#if defined SPARSE_FLAGS + +void DH_Subband::fillExtraData() +{ + BlobOStream& bos = createExtraBlob(); + + for (int stat = 0; stat < itsNrStations; stat ++) + itsFlags[stat].write(bos); +} + + +void DH_Subband::getExtraData() +{ + BlobIStream &bis = getExtraBlob(); + + for (int stat = 0; stat < itsNrStations; stat ++) + itsFlags[stat].read(bis); +} + +#endif + + +void DH_Subband::swapBytes() +{ + // only convert Samples; CEPframe converts Flags and Delays + dataConvert(LittleEndian, itsSamples, nrSamples()); +} +} // namespace CS1 } // namespace LOFAR diff --git a/Appl/CEP/CS1/CS1_Interface/src/Makefile.am b/Appl/CEP/CS1/CS1_Interface/src/Makefile.am index 9b45ee8b400c8eec39771e3e53b276e1099be201..c2384aaed4679e84707e07c6fc47fd19d8c13a5e 100644 --- a/Appl/CEP/CS1/CS1_Interface/src/Makefile.am +++ b/Appl/CEP/CS1/CS1_Interface/src/Makefile.am @@ -9,6 +9,7 @@ libcs1_interface_la_SOURCES = \ DH_Subband.cc \ DH_Visibilities.cc \ RSPTimeStamp.cc \ + SparseSet.cc \ Stub_BGL.cc \ Stub_BGL_RFI_Mitigation.cc \ Stub_BGL_Subband.cc \ diff --git a/Appl/CEP/CS1/CS1_Interface/src/SparseSet.cc b/Appl/CEP/CS1/CS1_Interface/src/SparseSet.cc new file mode 100644 index 0000000000000000000000000000000000000000..00b42312eec502fe49ccffdc218f97eda81747a3 --- /dev/null +++ b/Appl/CEP/CS1/CS1_Interface/src/SparseSet.cc @@ -0,0 +1,202 @@ +//# SparseSet.h: portable <bitset> adaptation +//# +//# Copyright (C) 2006 +//# ASTRON (Netherlands Foundation for Research in Astronomy) +//# P.O.Box 2, 7990 AA Dwingeloo, The Netherlands, seg@astron.nl +//# +//# 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 Free Software Foundation; either version 2 of the License, or +//# (at your option) any later version. +//# +//# This program is distributed in the hope that it will be useful, +//# but WITHOUT ANY WARRANTY; without even the implied warranty of +//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +//# GNU General Public License for more details. +//# +//# You should have received a copy of the GNU General Public License +//# along with this program; if not, write to the Free Software +//# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +//# +//# $Id$ + + +// #include "lofar_config.h" + +#include <CS1_Interface/SparseSet.h> + +#include <cstring> +#include <iostream> + + +namespace LOFAR +{ + +SparseSet &SparseSet::include(unsigned first, unsigned last) +{ + ++ last; + + if (ranges.size() == 0) { + struct range range = { first, last }; + ranges.push_back(range); + } else { + std::vector<range>::iterator last_involved = ranges.end(); + + while (-- last_involved >= ranges.begin() && last < last_involved->begin) + ; + + std::vector<range>::iterator first_involved = last_involved + 1; + + while (first_involved > ranges.begin() && first <= first_involved[-1].end) + -- first_involved; + + if (first_involved > last_involved) { + // insert new range + struct range range = { first, last }; + ranges.insert(first_involved, range); + } else { + // merge with existing range(s) + first_involved->begin = std::min(first, first_involved->begin); + first_involved->end = std::max(last, last_involved->end); + + ranges.erase(first_involved + 1, last_involved + 1); + } + } + + return *this; +} + + +unsigned SparseSet::count() const +{ + unsigned count = 0; + + for (std::vector<range>::const_iterator it = ranges.begin(); it != ranges.end(); it ++) + count += it->end - it->begin; + + return count; +} + + +bool SparseSet::test(unsigned index) const +{ + std::vector<range>::const_iterator it = ranges.begin(); + + while (it != ranges.end() && index >= it->begin) + it ++; + + return it != ranges.begin() && index < it[-1].end; +} + + +SparseSet SparseSet::operator | (const SparseSet &other) const +{ + if (ranges.size() == 0) { + return other; + } else if (other.ranges.size() == 0) { + return *this; + } else { + SparseSet union_set; + std::vector<range>::const_iterator range1 = ranges.begin(), range2 = other.ranges.begin(); + + while (range1 != ranges.end() && range2 != other.ranges.end()) { + if (range1->end < range2->begin) { + union_set.ranges.push_back(*range1 ++); + } else if (range2->end < range1->begin) { + union_set.ranges.push_back(*range2 ++); + } else { // there is overlap, or range1 and range2 are contiguous + struct range new_range; + new_range.begin = std::min(range1->begin, range2->begin); + + // check if subsequent ranges from set1 and set2 must be joined as well + while (1) { + if (range1 + 1 != ranges.end() && range1[1].begin <= range2->end) { + ++ range1; + } else if (range2 + 1 != other.ranges.end() && range2[1].begin <= range1->end) { + ++ range2; + } else { + break; + } + } + + new_range.end = std::max(range1->end, range2->end); + union_set.ranges.push_back(new_range); + ++ range1, ++ range2; + } + } + + union_set.ranges.insert(union_set.ranges.end(), range1, ranges.end()); + union_set.ranges.insert(union_set.ranges.end(), range2, other.ranges.end()); + return union_set; + } +} + + +SparseSet &SparseSet::operator -= (size_t count) +{ + std::vector<range>::iterator src = ranges.begin(), dst = src; + + while (src != ranges.end() && src->end <= count) + ++ src; + + if (src->begin < count) + src->begin = count; + + for (; src != ranges.end(); src ++, dst ++) { + dst->begin = src->begin - count; + dst->end = src->end - count; + } + + ranges.resize(ranges.size() + dst - src); + + return *this; +} + + +void SparseSet::write(BlobOStream &bos) const +{ + bos << (uint32) ranges.size(); + + for (std::vector<range>::const_iterator it = ranges.begin(); it != ranges.end(); it ++) + bos << (uint32) it->begin << (uint32) it->end; +} + + +void SparseSet::read(BlobIStream &bis) +{ + uint32 size, begin, end; + + bis >> size; + ranges.resize(size); + + for (std::vector<SparseSet::range>::iterator it = ranges.begin(); it != ranges.end(); it ++) { + bis >> begin >> end; + it->begin = begin; + it->end = end; + } +} + + + +#if 0 +SparseSet SparseSet::subset(unsigned first, unsigned last) const +{ + SparseSet subset; + size_t size = ranges.size(); + + if (size > 0) + for (const struct range *range = &ranges[0], *range_end = range + size; range < range_end; range ++) + if (first < range->end && last > range->begin) { + struct range new_range = { + std::max(first, range->begin) - first, + std::min(last, range->end) - first + }; + + subset.ranges.push_back(new_range); + } + + return subset; +} +#endif + +}