From f040b2f96fb6f3d95b73e9733cd3aa5e9044a140 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Andr=C3=A9=20Offringa?= <offringa@astron.nl>
Date: Sun, 9 Oct 2011 09:41:35 +0000
Subject: [PATCH] Bug 1491: allowing reading the mulitbeam raw files of the TKP

---
 .../include/AOFlagger/msio/rawdescfile.h      | 19 +++++
 .../include/AOFlagger/msio/rawreader.h        | 78 ++++++++++---------
 .../strategy/imagesets/rawdescimageset.h      |  8 +-
 3 files changed, 66 insertions(+), 39 deletions(-)

diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawdescfile.h b/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawdescfile.h
index d5b24b12992..fe1db726e2c 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawdescfile.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawdescfile.h
@@ -49,6 +49,12 @@ class RawDescFile
 		
 		unsigned TimestepsPerBlockCount() const { return _timestepsPerBlockCount; }
 		
+		unsigned BlockHeaderSize() const { return _blockHeaderSize; }
+		
+		unsigned BlockFooterSize() const { return _blockFooterSize; }
+		
+		unsigned SelectedBeam() const { return _selectedBeam; }
+		
 		const double TimeResolution() const { return _timeRes; }
 
 		double DisplayedTimeDuration() const { return _displayedTimeDuration; }
@@ -65,6 +71,11 @@ class RawDescFile
 		unsigned _channelsPerSubbandCount;
 		unsigned _timestepsPerBlockCount;
 		
+		unsigned _blockHeaderSize;
+		unsigned _blockFooterSize;
+		
+		unsigned _selectedBeam;
+		
 		double _timeRes;
 		double _displayedTimeDuration;
 		
@@ -85,6 +96,14 @@ class RawDescFile
 			std::getline(file, l);
 			_timestepsPerBlockCount = (unsigned) atol(l.c_str());
 			
+			std::getline(file, l);
+			_blockHeaderSize = (unsigned) atol(l.c_str());
+			std::getline(file, l);
+			_blockFooterSize = (unsigned) atol(l.c_str());
+			
+			std::getline(file, l);
+			_selectedBeam = (unsigned) atol(l.c_str());
+			
 			std::getline(file, l);
 			_timeRes = atof(l.c_str());
 			std::getline(file, l);
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawreader.h b/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawreader.h
index ec005c0a15c..2b35283b0ea 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawreader.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/msio/rawreader.h
@@ -31,11 +31,11 @@ class RawReader
 	public:
 		explicit RawReader(const std::string &filename) :
 		_blockHeaderSize(512),
-		_blockPostfixSize(8),
+		_blockFooterSize(8),
 		_beamCount(1),
 		_subbandCount(1),
 		_channelsPerSubbandCount(1),
-		_samplesPerBlockCount(12208),
+		_timestepsPerBlockCount(12208),
 		_filename(filename),
 		_useNetworkOrder(true),
 		_outputStream(0),
@@ -57,16 +57,16 @@ class RawReader
 			stream.seekg(0, std::ios_base::end);
 			std::streampos fileSize = stream.tellg();
 			unsigned long blockSize = BlockSize();
-			return (fileSize / blockSize) * _samplesPerBlockCount;
+			return (fileSize / blockSize) * _timestepsPerBlockCount;
 		}
 		
 		void Read(size_t startIndex, size_t endIndex, float *dest)
 		{
-			AOLogger::Debug << "Reading " << startIndex << " to " << endIndex << " (total: " << TimestepCount() << ")\n";
+			AOLogger::Debug << "Reading " << startIndex << " to " << endIndex << " (total: " << TimestepCount() << ") with " << _beamCount << " beams. \n";
 			
 			std::ifstream stream(_filename.c_str());
 			
-			size_t startBlock = startIndex / _samplesPerBlockCount;
+			size_t startBlock = startIndex / _timestepsPerBlockCount;
 			stream.seekg(startBlock * BlockSize(), std::ios_base::beg);
 			
 			const size_t rowsPerTimestep =  _beamCount * _subbandCount * _channelsPerSubbandCount;
@@ -75,22 +75,22 @@ class RawReader
 			RawBlock block(*this);
 			block.read(stream);
 			
-			size_t startTimestepInFirstBlock = startIndex - startBlock * _samplesPerBlockCount;
+			size_t startTimestepInFirstBlock = startIndex - startBlock * _timestepsPerBlockCount;
 			
-			size_t totalTimestepsInFirstBlock = (_samplesPerBlockCount - startTimestepInFirstBlock);
+			size_t totalTimestepsInFirstBlock = (_timestepsPerBlockCount - startTimestepInFirstBlock);
 			if(totalTimestepsInFirstBlock > endIndex - startIndex)
 				totalTimestepsInFirstBlock = endIndex - startIndex;
 			
-			memcpy(dest, block.SamplePtr(0, 0, 0, startTimestepInFirstBlock), bytesPerTimestep * totalTimestepsInFirstBlock);
+			memcpy(dest, block.SamplePtr(startTimestepInFirstBlock), bytesPerTimestep * totalTimestepsInFirstBlock);
 			
 			size_t samplesCopied = rowsPerTimestep * totalTimestepsInFirstBlock;
 			
-			size_t currentIndex = (startBlock + 1) * _samplesPerBlockCount;
+			size_t currentIndex = (startBlock + 1) * _timestepsPerBlockCount;
 			while(currentIndex < endIndex)
 			{
-				//AOLogger::Debug << currentIndex << '\n';
+				AOLogger::Debug << currentIndex << ", stream=" << stream.tellg() << "\n";
 				block.read(stream);
-				if(currentIndex + _samplesPerBlockCount > endIndex)
+				if(currentIndex + _timestepsPerBlockCount > endIndex)
 				{
 					// Whole block won't fit
 					memcpy(&dest[samplesCopied], block.SamplePtr(), bytesPerTimestep * (endIndex - currentIndex));
@@ -99,11 +99,12 @@ class RawReader
 				}
 				else {
 					// Block fits
-					memcpy(&dest[samplesCopied], block.SamplePtr(), bytesPerTimestep * _samplesPerBlockCount);
-					samplesCopied += rowsPerTimestep * _samplesPerBlockCount;
-					currentIndex += _samplesPerBlockCount;
+					memcpy(&dest[samplesCopied], block.SamplePtr(), bytesPerTimestep * _timestepsPerBlockCount);
+					samplesCopied += rowsPerTimestep * _timestepsPerBlockCount;
+					currentIndex += _timestepsPerBlockCount;
 				}
 			}
+			AOLogger::Debug << "Done reading.\n";
 		}
 		
 		void StartWrite()
@@ -116,9 +117,9 @@ class RawReader
 		{
 			if(_timeStepsInWriteBlock != 0)
 			{
-				for(size_t i=_timeStepsInWriteBlock;i<_samplesPerBlockCount;++i)
+				for(size_t i=_timeStepsInWriteBlock;i<_timestepsPerBlockCount;++i)
 				{
-					float *currentSamplePtr = _block.SamplePtr(0, 0, 0, i);
+					float *currentSamplePtr = _block.SamplePtr(i);
 					for(size_t j=0;j<_beamCount * _subbandCount * _channelsPerSubbandCount;++j)
 					{
 						currentSamplePtr[j] = 0.0;
@@ -132,11 +133,11 @@ class RawReader
 		
 		void Write(float *data, size_t timestepCount)
 		{
-			float *currentSamplePtr = _block.SamplePtr(0, 0, 0, _timeStepsInWriteBlock);
+			float *currentSamplePtr = _block.SamplePtr(_timeStepsInWriteBlock);
 			size_t writeCount;
 			
 			// Fill the current block
-			size_t writeTimesteps = _samplesPerBlockCount - _timeStepsInWriteBlock;
+			size_t writeTimesteps = _timestepsPerBlockCount - _timeStepsInWriteBlock;
 			if(writeTimesteps > timestepCount)
 				writeTimesteps = timestepCount;
 			writeCount = _beamCount * _subbandCount * _channelsPerSubbandCount * writeTimesteps;
@@ -146,24 +147,24 @@ class RawReader
 			_timeStepsInWriteBlock += writeTimesteps;
 			
 			// Has the block been filled?
-			if(_timeStepsInWriteBlock == _samplesPerBlockCount)
+			if(_timeStepsInWriteBlock == _timestepsPerBlockCount)
 			{
 				_block.write(*_outputStream);
 				_timeStepsInWriteBlock = 0;
 				
 				// Write whole blocks until no whole block of data is available
-				writeCount = _beamCount * _subbandCount * _channelsPerSubbandCount * _samplesPerBlockCount;
-				while(timestepCount >= _samplesPerBlockCount)
+				writeCount = _beamCount * _subbandCount * _channelsPerSubbandCount * _timestepsPerBlockCount;
+				while(timestepCount >= _timestepsPerBlockCount)
 				{
 					memcpy(_block.SamplePtr(), data, sizeof(float) * writeCount);
 					_block.write(*_outputStream);
-					timestepCount -= _samplesPerBlockCount;
-					data += _samplesPerBlockCount;
+					timestepCount -= _timestepsPerBlockCount;
+					data += _timestepsPerBlockCount;
 				}
 				
 				// Fill the last block as far as available
 				writeCount = _beamCount * _subbandCount * _channelsPerSubbandCount * timestepCount;
-				memcpy(_block.SamplePtr(0, 0, 0, _timeStepsInWriteBlock), data, sizeof(float) * writeCount);
+				memcpy(_block.SamplePtr(_timeStepsInWriteBlock), data, sizeof(float) * writeCount);
 				_timeStepsInWriteBlock += timestepCount;
 			}
 		}
@@ -172,20 +173,24 @@ class RawReader
 		
 		size_t BlockSize() const
 		{
-			return _blockHeaderSize + _blockPostfixSize +
-				_beamCount * _subbandCount * _channelsPerSubbandCount * _samplesPerBlockCount * sizeof(float);
+			return _blockHeaderSize + _blockFooterSize +
+				_beamCount * _subbandCount * _channelsPerSubbandCount * _timestepsPerBlockCount * sizeof(float);
 		}
 		
 		void SetSubbandCount(unsigned count) { _subbandCount = count; }
 		void SetChannelCount(unsigned count) { _channelsPerSubbandCount = count; }
+		void SetBeamCount(unsigned count) { _beamCount = count; }
+		void SetTimestepsPerBlockCount(unsigned count) { _timestepsPerBlockCount = count; }
+		void SetBlockHeaderSize(unsigned size) { _blockHeaderSize = size; }
+		void SetBlockFooterSize(unsigned size) { _blockFooterSize = size; }
 		
 	private:
 		unsigned _blockHeaderSize;
-		unsigned _blockPostfixSize;
+		unsigned _blockFooterSize;
 		unsigned _beamCount;
 		unsigned _subbandCount;
 		unsigned _channelsPerSubbandCount;
-		unsigned _samplesPerBlockCount;
+		unsigned _timestepsPerBlockCount;
 		const std::string _filename;
 		bool _useNetworkOrder;
 		
@@ -217,8 +222,8 @@ class RawReader
 				void initialize()
 				{
 					_header = new unsigned char[_reader._blockHeaderSize];
-					_data = new float[_reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._samplesPerBlockCount];
-					_postFix = new unsigned char[_reader._blockPostfixSize];
+					_data = new float[_reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._timestepsPerBlockCount];
+					_postFix = new unsigned char[_reader._blockFooterSize];
 				}
 				
 				void deinitialize()
@@ -230,10 +235,10 @@ class RawReader
 				
 				void read(std::istream &stream)
 				{
-					size_t length = _reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._samplesPerBlockCount;
+					size_t length = _reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._timestepsPerBlockCount;
 					stream.read(reinterpret_cast<char*>(_header), _reader._blockHeaderSize);
 					stream.read(reinterpret_cast<char*>(_data), length * sizeof(float));
-					stream.read(reinterpret_cast<char*>(_postFix), _reader._blockPostfixSize);
+					stream.read(reinterpret_cast<char*>(_postFix), _reader._blockFooterSize);
 					
 					if(_reader._useNetworkOrder)
 					{
@@ -246,7 +251,7 @@ class RawReader
 				
 				void write(std::ostream &stream)
 				{
-					size_t length = _reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._samplesPerBlockCount;
+					size_t length = _reader._beamCount * _reader._subbandCount * _reader._channelsPerSubbandCount * _reader._timestepsPerBlockCount;
 					if(_reader._useNetworkOrder)
 					{
 						for(size_t i=0;i<length;++i)
@@ -257,7 +262,7 @@ class RawReader
 					
 					stream.write(reinterpret_cast<char*>(_header), _reader._blockHeaderSize);
 					stream.write(reinterpret_cast<char*>(_data), length * sizeof(float));
-					stream.write(reinterpret_cast<char*>(_postFix), _reader._blockPostfixSize);
+					stream.write(reinterpret_cast<char*>(_postFix), _reader._blockFooterSize);
 					
 					if(_reader._useNetworkOrder)
 					{
@@ -273,12 +278,9 @@ class RawReader
 					return _data;
 				}
 				
-				float *SamplePtr(size_t beamIndex, size_t subbandIndex, size_t channelIndex, size_t sampleIndex)
+				float *SamplePtr(size_t sampleIndex)
 				{
 					size_t dataIndex =
-						channelIndex +
-						subbandIndex * _reader._channelsPerSubbandCount +
-						beamIndex * _reader._channelsPerSubbandCount * _reader._subbandCount +
 						sampleIndex * _reader._beamCount * _reader._channelsPerSubbandCount * _reader._subbandCount;
 					return &_data[dataIndex];
 				}
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/strategy/imagesets/rawdescimageset.h b/CEP/DP3/AOFlagger/include/AOFlagger/strategy/imagesets/rawdescimageset.h
index d7abc04297c..e8b68db1d81 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/strategy/imagesets/rawdescimageset.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/strategy/imagesets/rawdescimageset.h
@@ -80,6 +80,7 @@ namespace rfiStrategy {
 			{
 				AOLogger::Debug
 					<< "Opening rawdescfile, beams=" << _rawDescFile.BeamCount()
+					<< "(" << _rawDescFile.SelectedBeam() << ")"
 					<< ", subbands=" << _rawDescFile.SubbandCount()
 					<< ", channels=" << _rawDescFile.ChannelsPerSubbandCount()
 					<< ", timesteps=" << _rawDescFile.TimestepsPerBlockCount() << '\n';
@@ -91,6 +92,10 @@ namespace rfiStrategy {
 					_readers[i] = new RawReader(_rawDescFile.GetSet(i));
 					_readers[i]->SetSubbandCount(_rawDescFile.SubbandCount());
 					_readers[i]->SetChannelCount(_rawDescFile.ChannelsPerSubbandCount());
+					_readers[i]->SetBeamCount(_rawDescFile.BeamCount());
+					_readers[i]->SetTimestepsPerBlockCount(_rawDescFile.TimestepsPerBlockCount());
+					_readers[i]->SetBlockHeaderSize(_rawDescFile.BlockHeaderSize());
+					_readers[i]->SetBlockFooterSize(_rawDescFile.BlockFooterSize());
 					if(i == 0)
 						_totalTimesteps = _readers[i]->TimestepCount();
 					else if(_readers[i]->TimestepCount() < _totalTimesteps)
@@ -140,9 +145,10 @@ namespace rfiStrategy {
 				for(size_t setIndex=0;setIndex<_rawDescFile.GetCount();++setIndex)
 				{
 					_readers[setIndex]->Read(readStart, readStart + _imageWidth, data);
-					size_t pos = 0 * samplesPerTimestep; /*this selects beam zero for now*/
 					for(size_t x=0;x<_imageWidth;++x)
 					{
+						size_t pos = x * _rawDescFile.BeamCount() * _rawDescFile.ChannelsPerSubbandCount() * _rawDescFile.SubbandCount()
+							+ _rawDescFile.SelectedBeam() * _rawDescFile.ChannelsPerSubbandCount() * _rawDescFile.SubbandCount();
 						size_t y = setIndex * _rawDescFile.ChannelsPerSubbandCount() * _rawDescFile.SubbandCount();
 						for(size_t subbandIndex=0;subbandIndex < _rawDescFile.SubbandCount();++subbandIndex)
 						{
-- 
GitLab