Спектрограммы • dernasherbrezon

Введение

Мой проект r2cloud записывает сигнал пролетающего спутника на диск, а затем не спеша пытается извлечь из него данные и поделиться ими со всем миром. В большинстве случаев никаких данных получить не удаётся. Причины могут быть разными:

низкий пролёт спутника и, как следствие, слабый сигнал;
проблемы в приёмнике и/или антенне;
устаревшие параметры спутника (TLE);
спутник просто ничего не передаёт — либо нечего передавать, либо произошла ошибка, либо трансивер выключен.

Чтобы понять, в чём может быть проблема, удобно посмотреть на спектрограмму сигнала. Спектрограмма — это самый простой и при этом один из самых эффективных способов анализа.

Обычно спектрограмма представляет собой изображение. По оси X отложены частоты, по оси Y — время. Цветом кодируется мощность сигнала: синий или голубой — шум (сигнала нет), жёлтый — слабый сигнал, красный — сильный.

Такое представление позволяет быстро понять, был ли сигнал в определённое время на определённой частоте.

На изображении выше, например, хорошо видно:

три разных сигнала;
центральный сигнал — самый сильный;
все три сигнала смещаются справа налево. Это означает, что источник быстро движется, и частота меняется из-за эффекта Доплера;
сигналы появляются лишь периодически. Скорее всего, спутники передают телеметрию по таймеру.

Вычисление спектограммы

Чтобы построить такое изображение, сначала нужно определить, сколько герц будет приходиться на один пиксель:

int width = (int) (sampleRate / numHertzPerPixel);

Затем выполнить одномерное преобразование Фурье:

FloatFFT_1D fft = new FloatFFT_1D(width);
float[] complexBuf = new float[width * 2];
for (int i = 0; i < complexBuf.length; i += 2) {
	complexBuf[i] = source.readFloat();
	complexBuf[i + 1] = source.readFloat();
}
fft.complexForward(complexBuf);

Нормализовать результаты относительно ширины изображения:

float iNormalizationFactor = (float) 1 / width;
for (int i = 0, j = 0; i < complexBuf.length; i += 2, j++) {
	complexBuf[i] = complexBuf[i] * iNormalizationFactor;
	complexBuf[i + 1] = complexBuf[i + 1] * iNormalizationFactor;
}

Посчитать мощность:

for (int i = 0, j = 0; i < complexBuf.length; i += 2, j++) {
	tempResults[j] = (float) (10.0 * Math.log10((real * real) + (img * img) + 1e-20));
}

В итоге получается довольно длинное изображение. Чтобы этого избежать, можно принять один пиксель по высоте равным одной секунде и накапливать максимальное значение:

for (int j = 0; j < tempResults.length; j++) {
	tempResults[j] = Math.max(tempResults[j], (float) (10.0 * Math.log10((real * real) + (img * img) + 1e-20)));
}

Обычно середина изображения соответствует 0 Гц (baseband), поэтому две половинки нужно поменять местами:

for (int i = 0; i < length; i++) {
	float temp = tempResults[currentRow * width + i];
	tempResults[currentRow * width + i] = tempResults[currentRow * width + length + i];
	tempResults[currentRow * width + length + i] = temp;
}

После этого можно либо отнормировать результаты к минимальному и максимальному значениям, либо сразу создавать изображение:

BufferedImage image = new BufferedImage(width, height, BufferedImage.TYPE_INT_RGB);
for (int i = 0; i < height; i++) {
	for (int j = 0; j < width; j++) {
		image.setRGB(j, i, convertGrayToRgb(tempResults[width * i + j]));
	}
}

Алгоритм простой и вполне достаточный для анализа сигналов. Однако у него есть несколько недостатков, которые хотелось бы исправить:

PNG хранит 3 байта на пиксель, хотя достаточно одного. Из-за этого спектрограмма 10-минутного наблюдения занимает около 1,2 МБ.
Нет доступа к исходным данным. Поэтому нельзя отмасштабировать спектрограмму или корректно сравнить результаты с разных SDR-устройств.

Решение

Решение лежит на поверхности:

генерировать и хранить чёрно-белое изображение;
при показе конвертировать его в цветное.

Java отлично поддерживает чёрно-белые типы изображений, но есть важный нюанс. Значения пикселей нужно выставлять напрямую, в обход ColorModel. Иначе пиксели будут “улучшены” для более приятного восприятия человеком. В нашем случае спектрограмма — инструмент анализа данных, и у неё другая цель.

BufferedImage image = new BufferedImage(width, height, BufferedImage.TYPE_BYTE_GRAY);
WritableRaster r = image.getRaster();
int[] pixelArray = new int[1];
for (int i = 0; i < height; i++) {
	for (int j = 0; j < width; j++) {
		int pixel = (int) (tempResults[i * width + j] + 255);
		pixelArray[0] = (pixel & 0xFF);
		r.setPixel(j, height - i - 1, pixelArray);
	}
}

При показе на странице изображение дополнительно конвертируется. Делается это примерно так (псевдокод):

vm.$http.get(this.src, { baseURL: null, responseType: 'arraybuffer' }).then(function (response) {
  vm.canvas = document.getElementById('canvas');
  const ctx = vm.canvas.getContext('2d');
  const header = vm.parsePngHeader(response.data);
  const img = new Image();
  img.onload = () => {
    vm.canvas.width = img.width;
    vm.canvas.height = img.height;
    ctx.drawImage(img, 0, 0);
    const imageData = ctx.getImageData(0, 0, img.width, img.height);
    const data = imageData.data;
		  
    for (let j = 0; j < data.length; j += 4) {
      let rgb = vm.convertGrayToRgb(data[j] - 255);
      data[j] = rgb[0];
      data[j+1] = rgb[1];
      data[j+2] = rgb[2];
      data[j+3] = 0xFF;
    }
    ctx.putImageData(imageData, 0, 0);
  };
  img.src = URL.createObjectURL(new Blob([response.data], { type: 'image/png' }));
}).catch(function (error) {
  console.log(error)
})

Из интересных особенностей:

Чтобы получить доступ к чёрно-белым пикселям, изображение сначала нужно отрисовать.
Safari прогоняет пиксели через свою собственную ColorModel, поэтому они могут чуть-чуть отличаться от тех, что были записаны в Java. В других браузерах ситуация может быть немного лучше или хуже.

Приятный бонус отрисовки в браузере — возможность добавить интерактивность. Например, показывать время, частоту и мощность в текущем выделенном пикселе: