11#include "fitsbahtinovdetector.h"
12#include "fitsthresholddetector.h"
13#include "fitsgradientdetector.h"
14#include "fitscentroiddetector.h"
15#include "fitssepdetector.h"
19#include "kstarsdata.h"
23#include "skymapcomposite.h"
24#include "auxiliary/ksnotification.h"
25#include "auxiliary/robuststatistics.h"
28#include <QApplication>
30#include <QtConcurrent>
31#include <QImageReader>
33#include <QNetworkAccessManager>
35#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
40#if !defined(KSTARS_LITE) && defined(HAVE_LIBRAW)
41#include <libraw/libraw.h>
51#include <fits_debug.h>
53#define ZOOM_DEFAULT 100.0
56#define ZOOM_LOW_INCR 10
57#define ZOOM_HIGH_INCR 50
65const QStringList RAWFormats = {
"cr2",
"cr3",
"crw",
"nef",
"raf",
"dng",
"arw",
"orf" };
67bool FITSData::readableFilename(
const QString &filename)
70 QString extension = info.completeSuffix().toLower();
77FITSData::FITSData(FITSMode fitsMode): m_Mode(fitsMode)
81 qRegisterMetaType<FITSMode>(
"FITSMode");
83 debayerParams.method = DC1394_BAYER_METHOD_NEAREST;
84 debayerParams.filter = DC1394_COLOR_FILTER_RGGB;
85 debayerParams.offsetX = debayerParams.offsetY = 0;
88 m_CumulativeFrequency.resize(3);
89 m_HistogramBinWidth.resize(3);
90 m_HistogramFrequency.resize(3);
91 m_HistogramIntensity.resize(3);
102 qRegisterMetaType<FITSMode>(
"FITSMode");
104 debayerParams.method = DC1394_BAYER_METHOD_NEAREST;
105 debayerParams.filter = DC1394_COLOR_FILTER_RGGB;
106 debayerParams.offsetX = debayerParams.offsetY = 0;
110 this->m_Mode = other->m_Mode;
111 this->m_Statistics.channels = other->m_Statistics.channels;
112 memcpy(&m_Statistics, &(other->m_Statistics),
sizeof(m_Statistics));
113 m_ImageBuffer =
new uint8_t[m_Statistics.samples_per_channel * m_Statistics.channels * m_Statistics.bytesPerPixel];
114 memcpy(m_ImageBuffer, other->m_ImageBuffer,
115 m_Statistics.samples_per_channel * m_Statistics.channels * m_Statistics.bytesPerPixel);
133 if (m_WCSHandle !=
nullptr)
135 wcsvfree(&m_nwcs, &m_WCSHandle);
136 m_WCSHandle =
nullptr;
141 if (starCenters.
count() > 0)
142 qDeleteAll(starCenters);
145 if (m_SkyObjects.
count() > 0)
146 qDeleteAll(m_SkyObjects);
147 m_SkyObjects.
clear();
149 m_CatObjects.
clear();
153 fits_flush_file(fptr, &status);
154 fits_close_file(fptr, &status);
156 m_PackBuffer =
nullptr;
161void FITSData::loadCommon(
const QString &inFilename)
164 qDeleteAll(starCenters);
169 fits_flush_file(fptr, &status);
170 fits_close_file(fptr, &status);
172 m_PackBuffer =
nullptr;
176 m_Filename = inFilename;
179bool FITSData::loadFromBuffer(
const QByteArray &buffer)
183 return privateLoad(buffer);
188 loadCommon(inFilename);
190 m_Extension = info.completeSuffix().toLower();
191 qCDebug(KSTARS_FITS) <<
"Loading file " << m_Filename;
192#if QT_VERSION >= QT_VERSION_CHECK(6, 0, 0)
202QString fitsErrorToString(
int status)
204 char error_status[512] = {0};
205 fits_report_error(stderr, status);
206 fits_get_errstatus(status, error_status);
207 QString message = error_status;
212bool FITSData::privateLoad(
const QByteArray &buffer)
216 cacheEccentricity = -1;
219 return loadFITSImage(buffer);
221 return loadXISFImage(buffer);
223 return loadCanonicalImage(buffer);
224 else if (RAWFormats.
contains(m_Extension))
225 return loadRAWImage(buffer);
230bool FITSData::loadFITSImage(
const QByteArray &buffer,
const bool isCompressed)
232 int status = 0, anynull = 0;
235 m_HistogramConstructed =
false;
237 if (m_Extension.
contains(
".fz") || isCompressed)
246 m_compressedFilename = m_Filename;
251 rc = fp_unpack_file_to_fits(m_Filename.
toLocal8Bit().
data(), &fptr, fpvar) == 0;
254 m_Filename = uncompressedFile;
259 size_t m_PackBufferSize = 100000;
261 m_PackBuffer = (uint8_t *)malloc(m_PackBufferSize);
262 rc = fp_unpack_data_to_data(buffer.
data(), buffer.
size(), &m_PackBuffer, &m_PackBufferSize, fpvar) == 0;
266 void *data =
reinterpret_cast<void *
>(m_PackBuffer);
267 if (fits_open_memfile(&fptr, m_Filename.
toLocal8Bit().
data(), READONLY, &data, &m_PackBufferSize, 0,
270 m_LastError =
i18n(
"Error reading fits buffer: %1.", fitsErrorToString(status));
274 m_Statistics.size = m_PackBufferSize;
282 m_PackBuffer =
nullptr;
283 m_LastError =
i18n(
"Failed to unpack compressed fits");
284 qCCritical(KSTARS_FITS) << m_LastError;
288 m_isTemporary =
true;
289 m_isCompressed =
true;
290 m_Statistics.size = fptr->Fptr->logfilesize;
297 if (fits_open_diskfile(&fptr, m_Filename.
toLocal8Bit(), READONLY, &status))
299 m_LastError =
i18n(
"Error opening fits file %1 : %2", m_Filename, fitsErrorToString(status));
302 m_Statistics.size =
QFile(m_Filename).
size();
307 void *temp_buffer =
const_cast<void *
>(
reinterpret_cast<const void *
>(buffer.
data()));
308 size_t temp_size = buffer.
size();
309 if (fits_open_memfile(&fptr, m_Filename.
toLocal8Bit().
data(), READONLY,
310 &temp_buffer, &temp_size, 0,
nullptr, &status))
312 m_LastError =
i18n(
"Error reading fits buffer: %1", fitsErrorToString(status));
316 m_Statistics.size = temp_size;
319 if (fits_movabs_hdu(fptr, 1, IMAGE_HDU, &status))
323 m_PackBuffer =
nullptr;
324 m_LastError =
i18n(
"Could not locate image HDU: %1", fitsErrorToString(status));
327 if (fits_get_img_param(fptr, 3, &m_FITSBITPIX, &(m_Statistics.ndim), naxes, &status))
330 m_PackBuffer =
nullptr;
331 m_LastError =
i18n(
"FITS file open error (fits_get_img_param): %1", fitsErrorToString(status));
336 if ((fits_is_compressed_image(fptr, &status) || m_Statistics.ndim <= 0) && !isCompressed)
338 loadCommon(m_Filename);
339 qCDebug(KSTARS_FITS) <<
"Image is compressed. Reloading...";
340 return loadFITSImage(buffer,
true);
343 if (m_Statistics.ndim < 2)
345 m_LastError =
i18n(
"1D FITS images are not supported in KStars.");
346 qCCritical(KSTARS_FITS) << m_LastError;
348 m_PackBuffer =
nullptr;
352 switch (m_FITSBITPIX)
355 m_Statistics.dataType = TBYTE;
356 m_Statistics.bytesPerPixel =
sizeof(uint8_t);
360 m_FITSBITPIX = USHORT_IMG;
361 m_Statistics.dataType = TUSHORT;
362 m_Statistics.bytesPerPixel =
sizeof(int16_t);
365 m_Statistics.dataType = TUSHORT;
366 m_Statistics.bytesPerPixel =
sizeof(uint16_t);
370 m_FITSBITPIX = ULONG_IMG;
371 m_Statistics.dataType = TULONG;
372 m_Statistics.bytesPerPixel =
sizeof(int32_t);
375 m_Statistics.dataType = TULONG;
376 m_Statistics.bytesPerPixel =
sizeof(uint32_t);
379 m_Statistics.dataType = TFLOAT;
380 m_Statistics.bytesPerPixel =
sizeof(float);
383 m_Statistics.dataType = TLONGLONG;
384 m_Statistics.bytesPerPixel =
sizeof(int64_t);
387 m_Statistics.dataType = TDOUBLE;
388 m_Statistics.bytesPerPixel =
sizeof(double);
391 m_LastError =
i18n(
"Bit depth %1 is not supported.", m_FITSBITPIX);
392 qCCritical(KSTARS_FITS) << m_LastError;
396 if (m_Statistics.ndim < 3)
399 if (naxes[0] == 0 || naxes[1] == 0)
401 m_LastError =
i18n(
"Image has invalid dimensions %1x%2", naxes[0], naxes[1]);
402 qCCritical(KSTARS_FITS) << m_LastError;
404 m_PackBuffer =
nullptr;
408 m_Statistics.width = naxes[0];
409 m_Statistics.height = naxes[1];
410 m_Statistics.samples_per_channel = m_Statistics.width * m_Statistics.height;
411 roiCenter.
setX(m_Statistics.width / 2);
412 roiCenter.
setY(m_Statistics.height / 2);
413 if(m_Statistics.width % 2)
414 roiCenter.
setX(roiCenter.
x() + 1);
415 if(m_Statistics.height % 2)
416 roiCenter.
setY(roiCenter.
y() + 1);
420 m_Statistics.channels = naxes[2];
424 if ( (m_Mode != FITS_NORMAL && m_Mode != FITS_CALIBRATE) || !Options::auto3DCube())
425 m_Statistics.channels = 1;
427 m_ImageBufferSize = m_Statistics.samples_per_channel * m_Statistics.channels * m_Statistics.bytesPerPixel;
428 m_ImageBuffer =
new uint8_t[m_ImageBufferSize];
429 if (m_ImageBuffer ==
nullptr)
431 qCWarning(KSTARS_FITS) <<
"FITSData: Not enough memory for image_buffer channel. Requested: "
432 << m_ImageBufferSize <<
" bytes.";
435 m_PackBuffer =
nullptr;
442 long nelements = m_Statistics.samples_per_channel * m_Statistics.channels;
444 if (fits_read_img(fptr, m_Statistics.dataType, 1, nelements,
nullptr, m_ImageBuffer, &anynull, &status))
446 m_LastError =
i18n(
"Error reading image: %1", fitsErrorToString(status));
454 if (getRecordValue(
"DATE-OBS", value) && value.
isValid())
462 if (naxes[2] == 1 && m_Statistics.channels == 1 && Options::autoDebayer() && checkDebayer())
465 if (m_isTemporary && m_TemporaryDataFile.
open())
467 m_TemporaryDataFile.
write(buffer);
468 m_TemporaryDataFile.
close();
469 m_Filename = m_TemporaryDataFile.
fileName();
473 calculateStats(
false,
false);
476 calculateStats(
false,
false);
478 if (m_Mode == FITS_NORMAL || m_Mode == FITS_ALIGN)
481 starsSearched =
false;
486bool FITSData::loadXISFImage(
const QByteArray &buffer)
488 m_HistogramConstructed =
false;
494 LibXISF::XISFReader xisfReader;
501 LibXISF::ByteArray byteArray(buffer.
constData(), buffer.
size());
502 xisfReader.open(byteArray);
505 if (xisfReader.imagesCount() == 0)
507 m_LastError =
i18n(
"File contain no images");
511 const LibXISF::Image &image = xisfReader.getImage(0);
513 switch (image.sampleFormat())
515 case LibXISF::Image::UInt8:
516 m_Statistics.dataType = TBYTE;
517 m_Statistics.bytesPerPixel =
sizeof(LibXISF::UInt8);
518 m_FITSBITPIX = TBYTE;
520 case LibXISF::Image::UInt16:
521 m_Statistics.dataType = TUSHORT;
522 m_Statistics.bytesPerPixel =
sizeof(LibXISF::UInt16);
523 m_FITSBITPIX = TUSHORT;
525 case LibXISF::Image::UInt32:
526 m_Statistics.dataType = TULONG;
527 m_Statistics.bytesPerPixel =
sizeof(LibXISF::UInt32);
528 m_FITSBITPIX = TULONG;
530 case LibXISF::Image::Float32:
531 m_Statistics.dataType = TFLOAT;
532 m_Statistics.bytesPerPixel =
sizeof(LibXISF::Float32);
533 m_FITSBITPIX = TFLOAT;
536 m_LastError =
i18n(
"Sample format %1 is not supported.", LibXISF::Image::sampleFormatString(image.sampleFormat()).c_str());
537 qCCritical(KSTARS_FITS) << m_LastError;
541 m_Statistics.width = image.width();
542 m_Statistics.height = image.height();
543 m_Statistics.samples_per_channel = m_Statistics.width * m_Statistics.height;
544 m_Statistics.channels = image.channelCount();
545 m_Statistics.size = buffer.
size();
546 roiCenter.
setX(m_Statistics.width / 2);
547 roiCenter.
setY(m_Statistics.height / 2);
548 if(m_Statistics.width % 2)
549 roiCenter.
setX(roiCenter.
x() + 1);
550 if(m_Statistics.height % 2)
551 roiCenter.
setY(roiCenter.
y() + 1);
553 m_HeaderRecords.clear();
554 auto &fitsKeywords = image.fitsKeywords();
555 for(
auto &fitsKeyword : fitsKeywords)
562 m_ImageBufferSize = image.imageDataSize();
563 m_ImageBuffer =
new uint8_t[m_ImageBufferSize];
564 std::memcpy(m_ImageBuffer, image.imageData(), m_ImageBufferSize);
566 calculateStats(
false,
false);
569 catch (LibXISF::Error &error)
571 m_LastError =
i18n(
"XISF file open error: ") +
error.what();
582bool FITSData::saveXISFImage(
const QString &newFilename)
587 LibXISF::XISFWriter xisfWriter;
588 LibXISF::Image image;
589 image.setGeometry(m_Statistics.width, m_Statistics.height, m_Statistics.channels);
590 if (m_Statistics.channels > 1)
591 image.setColorSpace(LibXISF::Image::RGB);
593 switch (m_FITSBITPIX)
596 image.setSampleFormat(LibXISF::Image::UInt8);
599 image.setSampleFormat(LibXISF::Image::UInt16);
602 image.setSampleFormat(LibXISF::Image::UInt32);
605 image.setSampleFormat(LibXISF::Image::Float32);
608 m_LastError =
i18n(
"Bit depth %1 is not supported.", m_FITSBITPIX);
609 qCCritical(KSTARS_FITS) << m_LastError;
613 std::memcpy(image.imageData(), m_ImageBuffer, m_ImageBufferSize);
614 for (
auto &fitsKeyword : m_HeaderRecords)
615 image.addFITSKeyword({fitsKeyword.key.toUtf8().data(), fitsKeyword.value.toString().toUtf8().data(), fitsKeyword.comment.toUtf8().data()});
617 xisfWriter.writeImage(image);
619 m_Filename = newFilename;
621 catch (LibXISF::Error &err)
623 m_LastError =
i18n(
"Error saving XISF image") + err.what();
628 Q_UNUSED(newFilename)
633bool FITSData::loadCanonicalImage(
const QByteArray &buffer)
640 qCCritical(KSTARS_FITS) <<
"Failed to open image.";
644 m_Statistics.size = buffer.
size();
650 qCCritical(KSTARS_FITS) <<
"Failed to open image.";
662 m_FITSBITPIX = BYTE_IMG;
663 switch (m_FITSBITPIX)
666 m_Statistics.dataType = TBYTE;
667 m_Statistics.bytesPerPixel =
sizeof(uint8_t);
671 m_Statistics.dataType = TUSHORT;
672 m_Statistics.bytesPerPixel =
sizeof(int16_t);
675 m_Statistics.dataType = TUSHORT;
676 m_Statistics.bytesPerPixel =
sizeof(uint16_t);
680 m_Statistics.dataType = TULONG;
681 m_Statistics.bytesPerPixel =
sizeof(int32_t);
684 m_Statistics.dataType = TULONG;
685 m_Statistics.bytesPerPixel =
sizeof(uint32_t);
688 m_Statistics.dataType = TFLOAT;
689 m_Statistics.bytesPerPixel =
sizeof(float);
692 m_Statistics.dataType = TLONGLONG;
693 m_Statistics.bytesPerPixel =
sizeof(int64_t);
696 m_Statistics.dataType = TDOUBLE;
697 m_Statistics.bytesPerPixel =
sizeof(double);
700 m_LastError =
QString(
"Bit depth %1 is not supported.").
arg(m_FITSBITPIX);
701 qCCritical(KSTARS_FITS) << m_LastError;
705 m_Statistics.width =
static_cast<uint16_t
>(imageFromFile.
width());
706 m_Statistics.height =
static_cast<uint16_t
>(imageFromFile.
height());
708 m_Statistics.samples_per_channel = m_Statistics.width * m_Statistics.height;
710 m_ImageBufferSize = m_Statistics.samples_per_channel * m_Statistics.channels *
static_cast<uint16_t
>
711 (m_Statistics.bytesPerPixel);
712 m_ImageBuffer =
new uint8_t[m_ImageBufferSize];
713 if (m_ImageBuffer ==
nullptr)
715 m_LastError =
i18n(
"FITSData: Not enough memory for image_buffer channel. Requested: %1 bytes ", m_ImageBufferSize);
716 qCCritical(KSTARS_FITS) << m_LastError;
721 if (m_Statistics.channels == 1)
723 memcpy(m_ImageBuffer, imageFromFile.
bits(), m_ImageBufferSize);
728 auto debayered_buffer =
reinterpret_cast<uint8_t *
>(m_ImageBuffer);
729 auto * original_bayered_buffer =
reinterpret_cast<uint8_t *
>(imageFromFile.
bits());
732 uint8_t * rBuff = debayered_buffer;
733 uint8_t * gBuff = debayered_buffer + (m_Statistics.width * m_Statistics.height);
734 uint8_t * bBuff = debayered_buffer + (m_Statistics.width * m_Statistics.height * 2);
736 int imax = m_Statistics.samples_per_channel * 4 - 4;
737 for (
int i = 0; i <= imax; i += 4)
739 *rBuff++ = original_bayered_buffer[i + 2];
740 *gBuff++ = original_bayered_buffer[i + 1];
741 *bBuff++ = original_bayered_buffer[i + 0];
745 m_HeaderRecords.clear();
748 calculateStats(
false,
false);
753bool FITSData::loadRAWImage(
const QByteArray &buffer)
755#if !defined(KSTARS_LITE) && !defined(HAVE_LIBRAW)
756 m_LastError =
i18n(
"Unable to find dcraw and cjpeg. Please install the required tools to convert CR2/NEF to JPEG.");
770 m_LastError =
i18n(
"Cannot open file %1: %2", m_Filename, libraw_strerror(ret));
771 RawProcessor.recycle();
780 if ((ret = RawProcessor.open_buffer(
const_cast<void *
>(
reinterpret_cast<const void *
>(buffer.
data())),
781 buffer.
size())) != LIBRAW_SUCCESS)
783 m_LastError =
i18n(
"Cannot open buffer: %1", libraw_strerror(ret));
784 RawProcessor.recycle();
788 m_Statistics.size = buffer.
size();
792 if ((ret = RawProcessor.unpack()) != LIBRAW_SUCCESS)
794 m_LastError =
i18n(
"Cannot unpack_thumb: %1", libraw_strerror(ret));
795 RawProcessor.recycle();
799 if ((ret = RawProcessor.dcraw_process()) != LIBRAW_SUCCESS)
801 m_LastError =
i18n(
"Cannot dcraw_process: %1", libraw_strerror(ret));
802 RawProcessor.recycle();
806 libraw_processed_image_t *image = RawProcessor.dcraw_make_mem_image(&ret);
807 if (ret != LIBRAW_SUCCESS)
809 m_LastError =
i18n(
"Cannot load to memory: %1", libraw_strerror(ret));
810 RawProcessor.recycle();
814 RawProcessor.recycle();
816 m_Statistics.bytesPerPixel = image->bits / 8;
818 if (m_Statistics.bytesPerPixel == 1)
819 m_Statistics.dataType = TBYTE;
821 m_Statistics.dataType = TUSHORT;
822 m_Statistics.width = image->width;
823 m_Statistics.height = image->height;
824 m_Statistics.channels = image->colors;
825 m_Statistics.samples_per_channel = m_Statistics.width * m_Statistics.height;
827 m_ImageBufferSize = m_Statistics.samples_per_channel * m_Statistics.channels * m_Statistics.bytesPerPixel;
828 m_ImageBuffer =
new uint8_t[m_ImageBufferSize];
829 if (m_ImageBuffer ==
nullptr)
831 m_LastError =
i18n(
"FITSData: Not enough memory for image_buffer channel. Requested: %1 bytes ", m_ImageBufferSize);
832 qCCritical(KSTARS_FITS) << m_LastError;
833 libraw_dcraw_clear_mem(image);
838 auto destination_buffer =
reinterpret_cast<uint8_t *
>(m_ImageBuffer);
839 auto source_buffer =
reinterpret_cast<uint8_t *
>(image->data);
842 if (image->colors == 1)
844 memcpy(destination_buffer, source_buffer, m_ImageBufferSize);
849 uint8_t * rBuff = destination_buffer;
850 uint8_t * gBuff = destination_buffer + (m_Statistics.width * m_Statistics.height);
851 uint8_t * bBuff = destination_buffer + (m_Statistics.width * m_Statistics.height * 2);
853 int imax = m_Statistics.samples_per_channel * 3 - 3;
854 for (
int i = 0; i <= imax; i += 3)
856 *rBuff++ = source_buffer[i + 0];
857 *gBuff++ = source_buffer[i + 1];
858 *bBuff++ = source_buffer[i + 2];
861 libraw_dcraw_clear_mem(image);
863 m_HeaderRecords.clear();
866 calculateStats(
false,
false);
872bool FITSData::saveImage(
const QString &newFilename)
874 if (newFilename == m_Filename)
879 if (ext ==
"jpg" || ext ==
"png")
882 QImage fitsImage = FITSToImage(newFilename);
883 getMinMax(&min, &max);
889 else if (channels() == 1)
894 for (
int i = 0; i < 256; i++)
895 fitsImage.
setColor(i, qRgb(i, i, i));
902 double dataMin = m_Statistics.mean[0] - m_Statistics.stddev[0];
903 double dataMax = m_Statistics.mean[0] + m_Statistics.stddev[0] * 3;
905 double bscale = 255. / (dataMax - dataMin);
906 double bzero = (-dataMin) * (255. / (dataMax - dataMin));
909 switch (m_Statistics.dataType)
912 convertToQImage<uint8_t>(dataMin, dataMax, bscale, bzero, fitsImage);
916 convertToQImage<int16_t>(dataMin, dataMax, bscale, bzero, fitsImage);
920 convertToQImage<uint16_t>(dataMin, dataMax, bscale, bzero, fitsImage);
924 convertToQImage<int32_t>(dataMin, dataMax, bscale, bzero, fitsImage);
928 convertToQImage<uint32_t>(dataMin, dataMax, bscale, bzero, fitsImage);
932 convertToQImage<float>(dataMin, dataMax, bscale, bzero, fitsImage);
936 convertToQImage<int64_t>(dataMin, dataMax, bscale, bzero, fitsImage);
940 convertToQImage<double>(dataMin, dataMax, bscale, bzero, fitsImage);
949 m_Filename = newFilename;
950 qCInfo(KSTARS_FITS) <<
"Saved image file:" << m_Filename;
956 m_LastError =
i18n(
"Saving compressed files is not supported.");
968 fits_close_file(fptr, &status);
971 rotCounter = flipHCounter = flipVCounter = 0;
972 return saveXISFImage(newFilename);
1011 nelements = m_Statistics.samples_per_channel * m_Statistics.channels;
1014 if (fptr && fits_close_file(fptr, &status))
1016 m_LastError =
i18n(
"Failed to close file: %1", fitsErrorToString(status));
1021 if (fits_create_file(&new_fptr,
QString(
"!%1").arg(newFilename).toLocal8Bit(), &status))
1023 m_LastError =
i18n(
"Failed to create file: %1", fitsErrorToString(status));
1032 long naxis = m_Statistics.channels == 1 ? 2 : 3;
1033 long naxes[3] = {m_Statistics.width, m_Statistics.height, naxis};
1036 if (fits_create_img(fptr, m_FITSBITPIX, naxis, naxes, &status))
1038 m_LastError =
i18n(
"Failed to create image: %1", fitsErrorToString(status));
1045 if (fits_update_key(fptr, TDOUBLE,
"DATAMIN", &(m_Statistics.min),
"Minimum value", &status))
1047 m_LastError =
i18n(
"Failed to update key: %1", fitsErrorToString(status));
1052 if (fits_update_key(fptr, TDOUBLE,
"DATAMAX", &(m_Statistics.max),
"Maximum value", &status))
1054 m_LastError =
i18n(
"Failed to update key: %1", fitsErrorToString(status));
1059 fits_write_key(fptr, TDOUBLE,
"MIN1", &m_Statistics.min[0],
"Min Channel 1", &status);
1062 fits_write_key(fptr, TDOUBLE,
"MIN2", &m_Statistics.min[1],
"Min Channel 2", &status);
1063 fits_write_key(fptr, TDOUBLE,
"MIN3", &m_Statistics.min[2],
"Min Channel 3", &status);
1067 fits_write_key(fptr, TDOUBLE,
"MAX1", &m_Statistics.max[0],
"Max Channel 1", &status);
1070 fits_write_key(fptr, TDOUBLE,
"MAX2", &m_Statistics.max[1],
"Max Channel 2", &status);
1071 fits_write_key(fptr, TDOUBLE,
"MAX3", &m_Statistics.max[2],
"Max Channel 3", &status);
1075 if (m_Statistics.mean[0] > 0)
1077 fits_write_key(fptr, TDOUBLE,
"MEAN1", &m_Statistics.mean[0],
"Mean Channel 1", &status);
1080 fits_write_key(fptr, TDOUBLE,
"MEAN2", &m_Statistics.mean[1],
"Mean Channel 2", &status);
1081 fits_write_key(fptr, TDOUBLE,
"MEAN3", &m_Statistics.mean[2],
"Mean Channel 3", &status);
1086 if (m_Statistics.median[0] > 0)
1088 fits_write_key(fptr, TDOUBLE,
"MEDIAN1", &m_Statistics.median[0],
"Median Channel 1", &status);
1091 fits_write_key(fptr, TDOUBLE,
"MEDIAN2", &m_Statistics.median[1],
"Median Channel 2", &status);
1092 fits_write_key(fptr, TDOUBLE,
"MEDIAN3", &m_Statistics.median[2],
"Median Channel 3", &status);
1097 if (m_Statistics.stddev[0] > 0)
1099 fits_write_key(fptr, TDOUBLE,
"STDDEV1", &m_Statistics.stddev[0],
"Standard Deviation Channel 1", &status);
1102 fits_write_key(fptr, TDOUBLE,
"STDDEV2", &m_Statistics.stddev[1],
"Standard Deviation Channel 2", &status);
1103 fits_write_key(fptr, TDOUBLE,
"STDDEV3", &m_Statistics.stddev[2],
"Standard Deviation Channel 3", &status);
1108 for (
int i = 10; i < m_HeaderRecords.count(); i++)
1110 QString key = m_HeaderRecords[i].key;
1114 switch (value.
type())
1119 fits_write_key(fptr, TINT, key.
toLatin1().
constData(), &number, comment, &status);
1126 fits_write_key(fptr, TDOUBLE, key.
toLatin1().
constData(), &number, comment, &status);
1133 char valueBuffer[256] = {0};
1135 fits_write_key(fptr, TSTRING, key.
toLatin1().
constData(), valueBuffer, comment, &status);
1141 if (fits_write_date(fptr, &status))
1143 m_LastError =
i18n(
"Failed to update date: %1", fitsErrorToString(status));
1150 if (fits_write_history(fptr, history.
toLatin1(), &status))
1152 m_LastError =
i18n(
"Failed to update history: %1", fitsErrorToString(status));
1156 int rot = 0, mirror = 0;
1159 rot = (90 * rotCounter) % 360;
1163 if (flipHCounter % 2 != 0 || flipVCounter % 2 != 0)
1166 if ((rot != 0) || (mirror != 0))
1167 rotWCSFITS(rot, mirror);
1169 rotCounter = flipHCounter = flipVCounter = 0;
1172 if (fits_write_img(fptr, m_Statistics.dataType, 1, nelements, m_ImageBuffer, &status))
1174 m_LastError =
i18n(
"Failed to write image: %1", fitsErrorToString(status));
1178 m_Filename = newFilename;
1180 fits_flush_file(fptr, &status);
1182 qCInfo(KSTARS_FITS) <<
"Saved FITS file:" << m_Filename;
1187void FITSData::clearImageBuffers()
1189 delete[] m_ImageBuffer;
1190 m_ImageBuffer =
nullptr;
1191 if(m_ImageRoiBuffer !=
nullptr )
1193 delete[] m_ImageRoiBuffer;
1194 m_ImageRoiBuffer =
nullptr;
1200void FITSData::makeRoiBuffer(
QRect roi)
1205 uint32_t channelSize = roi.
height() * roi.
width();
1206 if(channelSize > m_Statistics.samples_per_channel || channelSize == 1)
1210 if(m_ImageRoiBuffer !=
nullptr )
1212 delete[] m_ImageRoiBuffer;
1213 m_ImageRoiBuffer =
nullptr;
1215 int xoffset = roi.
topLeft().
x() - 1;
1216 int yoffset = roi.
topLeft().
y() - 1;
1217 uint32_t bpp = m_Statistics.bytesPerPixel;
1218 m_ImageRoiBuffer =
new uint8_t[roi.
height()*roi.
width()*m_Statistics.channels * m_Statistics.bytesPerPixel]();
1219 for(
int n = 0 ; n < m_Statistics.channels ; n++)
1221 for(
int i = 0; i < roi.
height(); i++)
1223 size_t i1 = n * channelSize * bpp + i * roi.
width() * bpp;
1224 size_t i2 = n * m_Statistics.samples_per_channel * bpp + (yoffset + i) * width() * bpp + xoffset * bpp;
1225 memcpy(&m_ImageRoiBuffer[i1],
1231 memcpy(&m_ROIStatistics, &m_Statistics,
sizeof(FITSImage::Statistic));
1232 m_ROIStatistics.samples_per_channel = roi.
height() * roi.
width();
1233 m_ROIStatistics.width = roi.
width();
1234 m_ROIStatistics.height = roi.
height();
1235 calculateStats(
false,
true);
1238void FITSData::setupWCSParams()
1240 FITSImage::Solution solution;
1241 if (parseSolution(solution))
1243 const bool eastToTheRight = solution.parity == FITSImage::POSITIVE ? false :
true;
1244 updateWCSHeaderData(solution.orientation, solution.ra, solution.dec, solution.pixscale, eastToTheRight);
1247 bool validObservationDate =
false;
1249 if (getRecordValue(
"DATE-OBS", value))
1252 tsString = tsString.remove(
'\'').trimmed();
1260 validObservationDate =
true;
1266 if (!validObservationDate)
1271void FITSData::calculateStats(
bool refresh,
bool roi)
1276 calculateMinMax(refresh);
1277 calculateMedian(refresh);
1280 if (refresh ==
false && fptr)
1284 if (fits_read_key_dbl(fptr,
"MEAN1", &m_Statistics.mean[0],
nullptr, &status) == 0)
1287 fits_read_key_dbl(fptr,
"MEAN2", & m_Statistics.mean[1],
nullptr, &status);
1288 fits_read_key_dbl(fptr,
"MEAN3", &m_Statistics.mean[2],
nullptr, &status);
1291 if (fits_read_key_dbl(fptr,
"STDDEV1", &m_Statistics.stddev[0],
nullptr, &status) == 0)
1294 fits_read_key_dbl(fptr,
"STDDEV2", &m_Statistics.stddev[1],
nullptr, &status);
1295 fits_read_key_dbl(fptr,
"STDDEV3", &m_Statistics.stddev[2],
nullptr, &status);
1303 switch (m_Statistics.dataType)
1306 calculateStdDev<uint8_t>();
1310 calculateStdDev<int16_t>();
1314 calculateStdDev<uint16_t>();
1318 calculateStdDev<int32_t>();
1322 calculateStdDev<uint32_t>();
1326 calculateStdDev<float>();
1330 calculateStdDev<int64_t>();
1334 calculateStdDev<double>();
1342 m_Statistics.SNR = m_Statistics.mean[0] / m_Statistics.stddev[0];
1346 calculateMinMax(refresh, roi);
1347 calculateMedian(refresh, roi);
1349 switch (m_ROIStatistics.dataType)
1352 calculateStdDev<uint8_t>(roi);
1356 calculateStdDev<int16_t>(roi);
1360 calculateStdDev<uint16_t>(roi);
1364 calculateStdDev<int32_t>(roi);
1368 calculateStdDev<uint32_t>(roi);
1372 calculateStdDev<float>(roi);
1376 calculateStdDev<int64_t>(roi);
1380 calculateStdDev<double>(roi);
1389void FITSData::calculateMinMax(
bool refresh,
bool roi)
1399 if (fptr !=
nullptr && !refresh)
1404 if (fits_read_key_dbl(fptr,
"DATAMIN", &(m_Statistics.min[0]),
nullptr, &status) == 0)
1406 else if (fits_read_key_dbl(fptr,
"MIN1", &(m_Statistics.min[0]),
nullptr, &status) == 0)
1410 fits_read_key_dbl(fptr,
"MIN2", &m_Statistics.min[1],
nullptr, &status);
1411 fits_read_key_dbl(fptr,
"MIN3", &m_Statistics.min[2],
nullptr, &status);
1415 if (fits_read_key_dbl(fptr,
"DATAMAX", &(m_Statistics.max[0]),
nullptr, &status) == 0)
1417 else if (fits_read_key_dbl(fptr,
"MAX1", &(m_Statistics.max[0]),
nullptr, &status) == 0)
1421 fits_read_key_dbl(fptr,
"MAX2", &m_Statistics.max[1],
nullptr, &status);
1422 fits_read_key_dbl(fptr,
"MAX3", &m_Statistics.max[2],
nullptr, &status);
1425 if (nfound == 2 && !(m_Statistics.min[0] == 0 && m_Statistics.max[0] == 0))
1429 m_Statistics.min[0] = 1.0E30;
1430 m_Statistics.max[0] = -1.0E30;
1432 m_Statistics.min[1] = 1.0E30;
1433 m_Statistics.max[1] = -1.0E30;
1435 m_Statistics.min[2] = 1.0E30;
1436 m_Statistics.max[2] = -1.0E30;
1438 switch (m_Statistics.dataType)
1441 calculateMinMax<uint8_t>();
1445 calculateMinMax<int16_t>();
1449 calculateMinMax<uint16_t>();
1453 calculateMinMax<int32_t>();
1457 calculateMinMax<uint32_t>();
1461 calculateMinMax<float>();
1465 calculateMinMax<int64_t>();
1469 calculateMinMax<double>();
1478 m_ROIStatistics.min[0] = 1.0E30;
1479 m_ROIStatistics.max[0] = -1.0E30;
1481 m_ROIStatistics.min[1] = 1.0E30;
1482 m_ROIStatistics.max[1] = -1.0E30;
1484 m_ROIStatistics.min[2] = 1.0E30;
1485 m_ROIStatistics.max[2] = -1.0E30;
1487 switch (m_Statistics.dataType)
1490 calculateMinMax<uint8_t>(roi);
1494 calculateMinMax<int16_t>(roi);
1498 calculateMinMax<uint16_t>(roi);
1502 calculateMinMax<int32_t>(roi);
1506 calculateMinMax<uint32_t>(roi);
1510 calculateMinMax<float>(roi);
1514 calculateMinMax<int64_t>(roi);
1518 calculateMinMax<double>(roi);
1528void FITSData::calculateMedian(
bool refresh,
bool roi)
1538 if (fptr !=
nullptr && !refresh)
1541 if (fits_read_key_dbl(fptr,
"MEDIAN1", &m_Statistics.median[0],
nullptr, &status) == 0)
1545 fits_read_key_dbl(fptr,
"MEDIAN2", &m_Statistics.median[1],
nullptr, &status);
1546 fits_read_key_dbl(fptr,
"MEDIAN3", &m_Statistics.median[2],
nullptr, &status);
1552 m_Statistics.median[RED_CHANNEL] = 0;
1553 m_Statistics.median[GREEN_CHANNEL] = 0;
1554 m_Statistics.median[BLUE_CHANNEL] = 0;
1556 switch (m_Statistics.dataType)
1559 calculateMedian<uint8_t>();
1563 calculateMedian<int16_t>();
1567 calculateMedian<uint16_t>();
1571 calculateMedian<int32_t>();
1575 calculateMedian<uint32_t>();
1579 calculateMedian<float>();
1583 calculateMedian<int64_t>();
1587 calculateMedian<double>();
1596 m_ROIStatistics.median[RED_CHANNEL] = 0;
1597 m_ROIStatistics.median[GREEN_CHANNEL] = 0;
1598 m_ROIStatistics.median[BLUE_CHANNEL] = 0;
1600 switch (m_ROIStatistics.dataType)
1603 calculateMedian<uint8_t>(roi);
1607 calculateMedian<int16_t>(roi);
1611 calculateMedian<uint16_t>(roi);
1615 calculateMedian<int32_t>(roi);
1619 calculateMedian<uint32_t>(roi);
1623 calculateMedian<float>(roi);
1627 calculateMedian<int64_t>(roi);
1631 calculateMedian<double>(roi);
1641template <
typename T>
1642void FITSData::calculateMedian(
bool roi)
1644 auto * buffer =
reinterpret_cast<T *
>(roi ? m_ImageRoiBuffer : m_ImageBuffer);
1645 const uint32_t maxMedianSize = 500000;
1646 uint32_t medianSize = roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel;
1647 uint8_t downsample = 1;
1648 if (medianSize > maxMedianSize)
1650 downsample = (
static_cast<double>(medianSize) / maxMedianSize) + 0.999;
1651 medianSize /= downsample;
1657 std::vector<int32_t> samples;
1658 samples.reserve(medianSize);
1660 for (uint8_t n = 0; n < m_Statistics.channels; n++)
1662 auto *oneChannel = buffer + n * (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel);
1663 for (uint32_t upto = 0; upto < (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel);
1666 auto median = Mathematics::RobustStatistics::ComputeLocation(Mathematics::RobustStatistics::LOCATION_MEDIAN, samples);
1667 roi ? m_ROIStatistics.median[n] = median : m_Statistics.median[n] = median;
1671template <
typename T>
1672QPair<T, T> FITSData::getParitionMinMax(uint32_t start, uint32_t stride,
bool roi)
1674 auto * buffer =
reinterpret_cast<T *
>(roi ? m_ImageRoiBuffer : m_ImageBuffer);
1675 T min = std::numeric_limits<T>::max();
1676 T max = std::numeric_limits<T>::min();
1678 uint32_t
end = start + stride;
1680 for (uint32_t i = start; i <
end; i++)
1682 min = qMin(buffer[i], min);
1683 max = qMax(buffer[i], max);
1690 return qMakePair(min, max);
1693template <
typename T>
1694void FITSData::calculateMinMax(
bool roi)
1696 T min = std::numeric_limits<T>::max();
1697 T max = std::numeric_limits<T>::min();
1701 const uint8_t nThreads = 16;
1703 for (
int n = 0; n < m_Statistics.channels; n++)
1705 uint32_t cStart = n * (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel);
1708 uint32_t tStride = (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel) / nThreads;
1711 uint32_t fStride = tStride + ((roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel) -
1712 (tStride * nThreads));
1715 uint32_t tStart = cStart;
1720 for (
int i = 0; i < nThreads; i++)
1723#if QT_VERSION >= QT_VERSION_CHECK(6, 0, 0)
1724 futures.
append(
QtConcurrent::run(&FITSData::getParitionMinMax<T>,
this, tStart, (i == (nThreads - 1)) ? fStride : tStride,
1727 futures.
append(
QtConcurrent::run(
this, &FITSData::getParitionMinMax<T>, tStart, (i == (nThreads - 1)) ? fStride : tStride,
1734 for (
int i = 0; i < nThreads; i++)
1736 QPair<T, T> result = futures[i].result();
1737 min = qMin(result.first, min);
1738 max = qMax(result.second, max);
1743 m_Statistics.min[n] = min;
1744 m_Statistics.max[n] = max;
1748 m_ROIStatistics.min[n] = min;
1749 m_ROIStatistics.max[n] = max;
1762 SumData(
double s,
double sq,
int n) : sum(s), squaredSum(sq), numSamples(n) {}
1763 SumData() : sum(0), squaredSum(0), numSamples(0) {}
1766template <
typename T>
1767SumData getSumAndSquaredSum(uint32_t start, uint32_t stride, uint8_t *buff)
1769 auto * buffer =
reinterpret_cast<T *
>(buff);
1770 const uint32_t
end = start + stride;
1772 double squaredSum = 0;
1773 for (uint32_t i = start; i <
end; i++)
1775 double sample = buffer[i];
1777 squaredSum += sample * sample;
1779 const double numSamples =
end - start;
1780 return SumData(sum, squaredSum, numSamples);
1783template <
typename T>
1784void FITSData::calculateStdDev(
bool roi )
1787 const uint8_t nThreads = 16;
1789 for (
int n = 0; n < m_Statistics.channels; n++)
1791 uint32_t cStart = n * (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel);
1794 uint32_t tStride = (roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel) / nThreads;
1797 uint32_t fStride = tStride + ((roi ? m_ROIStatistics.samples_per_channel : m_Statistics.samples_per_channel) -
1798 (tStride * nThreads));
1801 uint32_t tStart = cStart;
1806 for (
int i = 0; i < nThreads; i++)
1809 uint8_t *buff = roi ? m_ImageRoiBuffer : m_ImageBuffer;
1811 (i == (nThreads - 1)) ? fStride : tStride, buff));
1816 double sum = 0, squared_sum = 0;
1817 double numSamples = 0;
1818 for (
int i = 0; i < nThreads; i++)
1820 const SumData result = futures[i].result();
1822 squared_sum += result.squaredSum;
1823 numSamples += result.numSamples;
1826 if (numSamples <= 0)
continue;
1827 const double mean = sum / numSamples;
1828 const double variance = squared_sum / numSamples - mean * mean;
1831 m_Statistics.mean[n] = mean;
1832 m_Statistics.stddev[n] = sqrt(variance);
1836 m_ROIStatistics.mean[n] = mean;
1837 m_ROIStatistics.stddev[n] = sqrt(variance);
1845 kernel.fill(0.0, size * size);
1847 double kernelSum = 0.0;
1848 int fOff = (size - 1) / 2;
1849 double normal = 1.0 / (2.0 * M_PI * sigma * sigma);
1850 for (
int y = -fOff; y <= fOff; y++)
1852 for (
int x = -fOff; x <= fOff; x++)
1854 double distance = ((y * y) + (x * x)) / (2.0 * sigma * sigma);
1855 int index = (y + fOff) * size + (x + fOff);
1856 kernel[index] = normal * qExp(-distance);
1857 kernelSum += kernel.at(index);
1860 for (
int y = 0; y < size; y++)
1862 for (
int x = 0; x < size; x++)
1864 int index = y * size + x;
1865 kernel[index] = kernel.at(index) * 1.0 / kernelSum;
1872template <
typename T>
1873void FITSData::convolutionFilter(
const QVector<double> &kernel,
int kernelSize)
1875 T * imagePtr =
reinterpret_cast<T *
>(m_ImageBuffer);
1881 int fOff = (kernelSize - 1) / 2;
1885 for (
int offsetY = 0; offsetY < m_Statistics.height; offsetY++)
1887 for (
int offsetX = 0; offsetX < m_Statistics.width; offsetX++)
1892 int byteOffset = offsetY * m_Statistics.width + offsetX;
1895 for (
int filterY = -fOff; filterY <= fOff; filterY++)
1897 for (
int filterX = -fOff; filterX <= fOff; filterX++)
1899 if ((offsetY + filterY) >= 0 && (offsetY + filterY) < m_Statistics.height
1900 && ((offsetX + filterX) >= 0 && (offsetX + filterX) < m_Statistics.width ))
1903 int calcOffset = byteOffset + filterX + filterY * m_Statistics.width;
1904 int index = (filterY + fOff) * kernelSize + (filterX + fOff);
1905 double kernelValue = kernel.
at(index);
1906 gt += (imagePtr[calcOffset]) * kernelValue;
1912 imagePtr[byteOffset] = gt;
1917template <
typename T>
1918void FITSData::gaussianBlur(
int kernelSize,
double sigma)
1921 if (kernelSize % 2 == 0)
1924 qCInfo(KSTARS_FITS) <<
"Warning, size must be an odd number, correcting size to " << kernelSize;
1932 QVector<double> gaussianKernel = createGaussianKernel(kernelSize, sigma);
1933 convolutionFilter<T>(gaussianKernel, kernelSize);
1936void FITSData::setMinMax(
double newMin,
double newMax, uint8_t channel)
1938 m_Statistics.min[channel] = newMin;
1939 m_Statistics.max[channel] = newMax;
1942bool FITSData::parseHeader()
1944 char * header =
nullptr;
1945 int status = 0, nkeys = 0;
1947 if (fits_hdr2str(fptr, 0,
nullptr, 0, &header, &nkeys, &status))
1949 fits_report_error(stderr, status);
1954 m_HeaderRecords.clear();
1957 for (
int i = 0; i < nkeys; i++)
1967 oneRecord.comment =
properties[0].mid(8).simplified();
1972 oneRecord.value =
properties[1].simplified();
1974 oneRecord.comment =
properties[2].simplified();
1981 oneRecord.value.toInt(&ok);
1987 oneRecord.value.toDouble(&ok);
1993 m_HeaderRecords.
append(oneRecord);
2001bool FITSData::getRecordValue(
const QString &key,
QVariant &value)
const
2003 auto result = std::find_if(m_HeaderRecords.begin(), m_HeaderRecords.end(), [&key](
const Record & oneRecord)
2005 return (oneRecord.key == key && oneRecord.value.isValid());
2008 if (result != m_HeaderRecords.end())
2010 value = (*result).
value;
2018 auto result = std::find_if(m_HeaderRecords.begin(), m_HeaderRecords.end(), [&key](
const Record & oneRecord)
2020 return (oneRecord.key == key && oneRecord.value.isValid());
2023 if (result != m_HeaderRecords.end())
2025 (*result).value = value;
2031#if QT_VERSION >= QT_VERSION_CHECK(6, 0, 0)
2032 m_HeaderRecords.
insert(std::max(0LL, m_HeaderRecords.size() - 1), record);
2034 m_HeaderRecords.
insert(std::max(0, m_HeaderRecords.size() - 1), record);
2039bool FITSData::parseSolution(FITSImage::Solution &solution)
const
2042 bool raOK =
false, deOK =
false, coordOK =
false, scaleOK =
false;
2046 solution.fieldWidth = solution.fieldHeight = solution.pixscale = solution.ra = solution.dec = -1;
2049 if (getRecordValue(
"OBJCTRA", value))
2052 solution.ra = angleValue.
Degrees();
2055 else if (getRecordValue(
"RA", value))
2057 solution.ra = value.
toDouble(&raOK);
2061 if (getRecordValue(
"OBJCTDEC", value))
2064 solution.dec = angleValue.
Degrees();
2067 else if (getRecordValue(
"DEC", value))
2069 solution.dec = value.
toDouble(&deOK);
2072 coordOK = raOK && deOK;
2076 if (getRecordValue(
"SCALE", value))
2081 double focal_length = -1;
2082 if (getRecordValue(
"FOCALLEN", value))
2087 double pixsize1 = -1, pixsize2 = -1;
2089 if (getRecordValue(
"PIXSIZE1", value))
2093 else if (getRecordValue(
"XPIXSZ", value))
2098 if (getRecordValue(
"PIXSIZE2", value))
2102 else if (getRecordValue(
"YPIXSZ", value))
2107 int binx = 1, biny = 1;
2109 if (getRecordValue(
"XBINNING", value))
2114 if (getRecordValue(
"YBINNING", value))
2119 if (pixsize1 > 0 && pixsize2 > 0)
2125 solution.pixscale = scale;
2127 solution.fieldWidth = (m_Statistics.width * scale) / 60.0;
2129 solution.fieldHeight = (m_Statistics.height * scale * (pixsize2 / pixsize1)) / 60.0;
2131 else if (focal_length > 0)
2134 solution.fieldWidth = ((206264.8062470963552 * m_Statistics.width * (pixsize1 / 1000.0)) / (focal_length * binx)) / 60.0;
2136 solution.fieldHeight = ((206264.8062470963552 * m_Statistics.height * (pixsize2 / 1000.0)) / (focal_length * biny)) / 60.0;
2138 solution.pixscale = (206264.8062470963552 * (pixsize1 / 1000.0)) / (focal_length * binx);
2144 return (coordOK || scaleOK);
2152 starAlgorithm = algorithm;
2153 qDeleteAll(starCenters);
2154 starCenters.
clear();
2155 starsSearched =
true;
2161 m_StarDetector.
reset(
new FITSSEPDetector(
this));
2162 m_StarDetector->setSettings(m_SourceExtractorSettings);
2163 if (m_Mode == FITS_NORMAL && trackingBox.
isNull())
2165 if (Options::quickHFR())
2168 const int w = getStatistics().width;
2169 const int h = getStatistics().height;
2170 QRect middle(
static_cast<int>(w * 0.25),
static_cast<int>(h * 0.25), w / 2, h / 2);
2171 m_StarFindFuture = m_StarDetector->findSources(middle);
2172 return m_StarFindFuture;
2175 m_StarFindFuture = m_StarDetector->findSources(trackingBox);
2176 return m_StarFindFuture;
2179 case ALGORITHM_GRADIENT:
2182 m_StarDetector.
reset(
new FITSGradientDetector(
this));
2183 m_StarDetector->setSettings(m_SourceExtractorSettings);
2184 m_StarFindFuture = m_StarDetector->findSources(trackingBox);
2185 return m_StarFindFuture;
2188 case ALGORITHM_CENTROID:
2191 m_StarDetector.
reset(
new FITSCentroidDetector(
this));
2192 m_StarDetector->setSettings(m_SourceExtractorSettings);
2194 if (!isHistogramConstructed())
2195 constructHistogram();
2196 m_StarDetector->configure(
"JMINDEX", m_JMIndex);
2197 m_StarFindFuture = m_StarDetector->findSources(trackingBox);
2198 return m_StarFindFuture;
2202 m_StarDetector.
reset(
new FITSCentroidDetector(
this));
2203 m_StarFindFuture = starDetector->findSources(trackingBox);
2204 return m_StarFindFuture;
2208 case ALGORITHM_THRESHOLD:
2210 m_StarDetector.
reset(
new FITSThresholdDetector(
this));
2211 m_StarDetector->setSettings(m_SourceExtractorSettings);
2212 m_StarDetector->configure(
"THRESHOLD_PERCENTAGE", Options::focusThreshold());
2213 m_StarFindFuture = m_StarDetector->findSources(trackingBox);
2214 return m_StarFindFuture;
2217 case ALGORITHM_BAHTINOV:
2219 m_StarDetector.
reset(
new FITSBahtinovDetector(
this));
2220 m_StarDetector->setSettings(m_SourceExtractorSettings);
2221 m_StarDetector->configure(
"NUMBER_OF_AVERAGE_ROWS", Options::focusMultiRowAverage());
2222 m_StarFindFuture = m_StarDetector->findSources(trackingBox);
2223 return m_StarFindFuture;
2230 if (mask.
isNull() ==
false)
2232 starCenters.
erase(std::remove_if(starCenters.
begin(), starCenters.
end(),
2235 return (mask->isVisible(edge->x, edge->y) == false);
2236 }), starCenters.
end());
2239 return starCenters.
count();
2243double FITSData::getHFR(HFRType type)
2245 if (starCenters.
empty())
2248 if (cacheHFR >= 0 && cacheHFRType == type)
2251 m_SelectedHFRStar.invalidate();
2253 if (type == HFR_MAX)
2258 for (
int i = 0; i < starCenters.
count(); i++)
2260 if (starCenters[i]->HFR > maxVal)
2263 maxVal = starCenters[i]->HFR;
2267 m_SelectedHFRStar = *starCenters[maxIndex];
2268 cacheHFR = starCenters[maxIndex]->HFR;
2269 cacheHFRType =
type;
2272 else if (type == HFR_HIGH)
2275 int minX = width() / 10;
2276 int minY = height() / 10;
2277 int maxX = width() - minX;
2278 int maxY = height() - minY;
2279 starCenters.
erase(std::remove_if(starCenters.
begin(), starCenters.
end(), [minX, minY, maxX, maxY](Edge * oneStar)
2281 return (oneStar->x < minX || oneStar->x > maxX || oneStar->y < minY || oneStar->y > maxY);
2282 }), starCenters.
end());
2284 if (starCenters.
empty())
2287 m_SelectedHFRStar = *starCenters[
static_cast<int>(starCenters.
size() * 0.05)];
2288 cacheHFR = m_SelectedHFRStar.HFR;
2289 cacheHFRType =
type;
2292 else if (type == HFR_MEDIAN)
2294 std::nth_element(starCenters.
begin(), starCenters.
begin() + starCenters.
size() / 2, starCenters.
end());
2295 m_SelectedHFRStar = *starCenters[starCenters.
size() / 2];
2297 cacheHFR = m_SelectedHFRStar.HFR;
2298 cacheHFRType =
type;
2306 const int saturationValue = m_Statistics.dataType == TBYTE ? 250 : 50000;
2307 int numSaturated = 0;
2308 for (
auto center : starCenters)
2309 if (
center->val > saturationValue)
2311 const bool removeSaturatedStars = starCenters.size() - numSaturated > 20;
2312 if (removeSaturatedStars && numSaturated > 0)
2313 qCDebug(KSTARS_FITS) <<
"Removing " << numSaturated <<
" stars from HFR calculation";
2315 std::vector<double> HFRs;
2317 for (
auto center : starCenters)
2319 if (removeSaturatedStars &&
center->val > saturationValue)
continue;
2321 if (type == HFR_AVERAGE)
2322 HFRs.push_back(
center->HFR);
2328 if (m_SkyBackground.mean <= 0.0 ||
center->val < m_SkyBackground.mean)
2330 HFRs.push_back(
center->HFR);
2331 qCDebug(KSTARS_FITS) <<
"HFR Adj, sky background " << m_SkyBackground.mean <<
" star peak " <<
center->val <<
2336 const double a_div_b =
center->val / m_SkyBackground.mean;
2337 const double factor = erf(sqrt(log(a_div_b))) / (1 - (1 / a_div_b));
2338 HFRs.push_back(
center->HFR * factor);
2339 qCDebug(KSTARS_FITS) <<
"HFR Adj, brightness adjusted from " <<
center->HFR <<
" to " <<
center->HFR * factor;
2344 auto m = Mathematics::RobustStatistics::ComputeLocation(Mathematics::RobustStatistics::LOCATION_SIGMACLIPPING,
2348 cacheHFRType = HFR_AVERAGE;
2352double FITSData::getHFR(
int x,
int y,
double scale)
2354 if (starCenters.empty())
2357 for (
int i = 0; i < starCenters.count(); i++)
2359 const int maxDist = std::max(2,
static_cast<int>(0.5 + 2 * starCenters[i]->width / scale));
2360 const int dx = std::fabs(starCenters[i]->x - x);
2361 const int dy = std::fabs(starCenters[i]->y - y);
2362 if (dx <= maxDist && dy <= maxDist)
2364 m_SelectedHFRStar = *starCenters[i];
2365 return starCenters[i]->HFR;
2372double FITSData::getEccentricity()
2374 if (starCenters.empty())
2376 if (cacheEccentricity >= 0)
2377 return cacheEccentricity;
2378 std::vector<float> eccs;
2379 for (
const auto &s : starCenters)
2380 eccs.push_back(s->ellipticity);
2381 int middle = eccs.size() / 2;
2382 std::nth_element(eccs.begin(), eccs.begin() + middle, eccs.end());
2383 float medianEllipticity = eccs[middle];
2389 const float eccentricity = sqrt(medianEllipticity * (2 - medianEllipticity));
2390 cacheEccentricity = eccentricity;
2391 return eccentricity;
2396 if (type == FITS_NONE)
2409 case FITS_AUTO_STRETCH:
2411 for (
int i = 0; i < 3; i++)
2413 dataMin[i] = m_Statistics.mean[i] - m_Statistics.stddev[i];
2414 dataMax[i] = m_Statistics.mean[i] + m_Statistics.stddev[i] * 3;
2419 case FITS_HIGH_CONTRAST:
2421 for (
int i = 0; i < 3; i++)
2423 dataMin[i] = m_Statistics.mean[i] + m_Statistics.stddev[i];
2424 dataMax[i] = m_Statistics.mean[i] + m_Statistics.stddev[i] * 3;
2429 case FITS_HIGH_PASS:
2431 for (
int i = 0; i < 3; i++)
2433 dataMin[i] = m_Statistics.mean[i];
2442 switch (m_Statistics.dataType)
2446 for (
int i = 0; i < 3; i++)
2448 dataMin[i] = dataMin[i] < 0 ? 0 : dataMin[i];
2449 dataMax[i] = dataMax[i] > UINT8_MAX ? UINT8_MAX : dataMax[i];
2451 applyFilter<uint8_t>(type, image, &dataMin, &dataMax);
2457 for (
int i = 0; i < 3; i++)
2459 dataMin[i] = dataMin[i] < INT16_MIN ? INT16_MIN : dataMin[i];
2460 dataMax[i] = dataMax[i] > INT16_MAX ? INT16_MAX : dataMax[i];
2462 applyFilter<uint16_t>(type, image, &dataMin, &dataMax);
2469 for (
int i = 0; i < 3; i++)
2471 dataMin[i] = dataMin[i] < 0 ? 0 : dataMin[i];
2472 dataMax[i] = dataMax[i] > UINT16_MAX ? UINT16_MAX : dataMax[i];
2474 applyFilter<uint16_t>(type, image, &dataMin, &dataMax);
2480 for (
int i = 0; i < 3; i++)
2482 dataMin[i] = dataMin[i] < INT_MIN ? INT_MIN : dataMin[i];
2483 dataMax[i] = dataMax[i] > INT_MAX ? INT_MAX : dataMax[i];
2485 applyFilter<uint16_t>(type, image, &dataMin, &dataMax);
2491 for (
int i = 0; i < 3; i++)
2493 dataMin[i] = dataMin[i] < 0 ? 0 : dataMin[i];
2494 dataMax[i] = dataMax[i] > UINT_MAX ? UINT_MAX : dataMax[i];
2496 applyFilter<uint16_t>(type, image, &dataMin, &dataMax);
2502 for (
int i = 0; i < 3; i++)
2504 dataMin[i] = dataMin[i] < FLT_MIN ? FLT_MIN : dataMin[i];
2505 dataMax[i] = dataMax[i] > FLT_MAX ? FLT_MAX : dataMax[i];
2507 applyFilter<float>(type, image, &dataMin, &dataMax);
2513 for (
int i = 0; i < 3; i++)
2515 dataMin[i] = dataMin[i] < LLONG_MIN ? LLONG_MIN : dataMin[i];
2516 dataMax[i] = dataMax[i] > LLONG_MAX ? LLONG_MAX : dataMax[i];
2519 applyFilter<long>(type, image, &dataMin, &dataMax);
2525 for (
int i = 0; i < 3; i++)
2527 dataMin[i] = dataMin[i] < DBL_MIN ? DBL_MIN : dataMin[i];
2528 dataMax[i] = dataMax[i] > DBL_MAX ? DBL_MAX : dataMax[i];
2530 applyFilter<double>(type, image, &dataMin, &dataMax);
2547template <
typename T>
2550 bool calcStats =
false;
2551 T * image =
nullptr;
2554 image =
reinterpret_cast<T *
>(targetImage);
2557 image =
reinterpret_cast<T *
>(m_ImageBuffer);
2562 for (
int i = 0; i < 3; i++)
2564 min[i] = (*targetMin)[i] < std::numeric_limits<T>::min() ? std::numeric_limits<T>::min() : (*targetMin)[i];
2565 max[i] = (*targetMax)[i] > std::numeric_limits<T>::max() ? std::numeric_limits<T>::max() : (*targetMax)[i];
2570 const uint8_t nThreads = 16;
2572 uint32_t width = m_Statistics.width;
2573 uint32_t height = m_Statistics.height;
2581 case FITS_AUTO_STRETCH:
2582 case FITS_HIGH_CONTRAST:
2585 case FITS_HIGH_PASS:
2591 if (type == FITS_LOG)
2593 for (
int i = 0; i < 3; i++)
2594 coeff[i] = max[i] / std::log(1 + max[i]);
2596 else if (type == FITS_SQRT)
2598 for (
int i = 0; i < 3; i++)
2599 coeff[i] = max[i] / sqrt(max[i]);
2602 for (
int n = 0; n < m_Statistics.channels; n++)
2604 if (type == FITS_HIGH_PASS)
2605 min[n] = m_Statistics.mean[n];
2607 uint32_t cStart = n * m_Statistics.samples_per_channel;
2610 uint32_t tStride = m_Statistics.samples_per_channel / nThreads;
2613 uint32_t fStride = tStride + (m_Statistics.samples_per_channel - (tStride * nThreads));
2615 T * runningBuffer = image + cStart;
2617 if (type == FITS_LOG)
2619 for (
int i = 0; i < nThreads; i++)
2622 futures.
append(
QtConcurrent::map(runningBuffer, (runningBuffer + ((i == (nThreads - 1)) ? fStride : tStride)), [min, max,
2625 a = qBound(min[n],
static_cast<T
>(round(coeff[n] * std::log(1 + qBound(min[n], a, max[n])))), max[n]);
2628 runningBuffer += tStride;
2631 else if (type == FITS_SQRT)
2633 for (
int i = 0; i < nThreads; i++)
2636 futures.
append(
QtConcurrent::map(runningBuffer, (runningBuffer + ((i == (nThreads - 1)) ? fStride : tStride)), [min, max,
2639 a = qBound(min[n],
static_cast<T
>(round(coeff[n] * a)), max[n]);
2643 runningBuffer += tStride;
2647 for (
int i = 0; i < nThreads; i++)
2650 futures.
append(
QtConcurrent::map(runningBuffer, (runningBuffer + ((i == (nThreads - 1)) ? fStride : tStride)), [min, max,
2653 a = qBound(min[n], a, max[n]);
2655 runningBuffer += tStride;
2660 for (
int i = 0; i < nThreads * m_Statistics.channels; i++)
2661 futures[i].waitForFinished();
2665 for (
int i = 0; i < 3; i++)
2667 m_Statistics.min[i] = min[i];
2668 m_Statistics.max[i] = max[i];
2670 calculateStdDev<T>();
2678 if (!isHistogramConstructed())
2679 constructHistogram();
2683 double coeff = 255.0 / (height * width);
2687 for (
int i = 0; i < m_Statistics.channels; i++)
2689 uint32_t offset = i * m_Statistics.samples_per_channel;
2690 for (uint32_t j = 0; j < height; j++)
2692 row = offset + j * width;
2693 for (uint32_t k = 0; k < width; k++)
2696 bufferVal = (image[index] - min[i]) / m_HistogramBinWidth[i];
2698 if (bufferVal >= m_CumulativeFrequency[i].size())
2699 bufferVal = m_CumulativeFrequency[i].size() - 1;
2701 image[index] = qBound(min[i],
static_cast<T
>(round(coeff * m_CumulativeFrequency[i][bufferVal])), max[i]);
2708 calculateStats(
true,
false);
2714 uint8_t BBP = m_Statistics.bytesPerPixel;
2715 auto * extension =
new T[(width + 2) * (height + 2)];
2720 for (uint32_t ch = 0; ch < m_Statistics.channels; ch++)
2722 uint32_t offset = ch * m_Statistics.samples_per_channel;
2723 uint32_t N = width, M = height;
2725 for (uint32_t i = 0; i < M; ++i)
2727 memcpy(extension + (N + 2) * (i + 1) + 1, image + (N * i) + offset, N * BBP);
2728 extension[(N + 2) * (i + 1)] = image[N * i + offset];
2729 extension[(N + 2) * (i + 2) - 1] = image[N * (i + 1) - 1 + offset];
2732 memcpy(extension, extension + N + 2, (N + 2) * BBP);
2734 memcpy(extension + (N + 2) * (M + 1), extension + (N + 2) * M, (N + 2) * BBP);
2741 for (uint32_t m = 1; m < M - 1; ++m)
2742 for (uint32_t n = 1; n < N - 1; ++n)
2748 memset(&window[0], 0, 9 *
sizeof(
float));
2749 for (uint32_t j = m - 1; j < m + 2; ++j)
2750 for (uint32_t i = n - 1; i < n + 2; ++i)
2751 window[k++] = extension[j * N + i];
2753 for (uint32_t j = 0; j < 5; ++j)
2757 for (uint32_t l = j + 1; l < 9; ++l)
2758 if (window[l] < window[mine])
2761 const float temp =
window[j];
2766 image[(m - 1) * (N - 2) + n - 1 + offset] =
window[4];
2774 calculateStdDev<T>();
2779 gaussianBlur<T>(Options::focusGaussianKernelSize(), Options::focusGaussianSigma());
2781 calculateStats(
true,
false);
2784 case FITS_ROTATE_CW:
2789 case FITS_ROTATE_CCW:
2794 case FITS_MOUNT_FLIP_H:
2799 case FITS_MOUNT_FLIP_V:
2812 for (
int i = 0; i < starCenters.count(); i++)
2814 int x =
static_cast<int>(starCenters[i]->x);
2815 int y =
static_cast<int>(starCenters[i]->y);
2818 starCentersInSubFrame.
append(starCenters[i]);
2821 return starCentersInSubFrame;
2824bool FITSData::loadWCS()
2826#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
2828 if (m_WCSState == Success)
2830 qCWarning(KSTARS_FITS) <<
"WCS data already loaded";
2834 if (m_WCSHandle !=
nullptr)
2836 wcsvfree(&m_nwcs, &m_WCSHandle);
2838 m_WCSHandle =
nullptr;
2841 qCDebug(KSTARS_FITS) <<
"Started WCS Data Processing...";
2845 int nkeyrec = 0, nreject = 0;
2848 char *header =
nullptr;
2849 if (fits_hdr2str(fptr, 1,
nullptr, 0, &header, &nkeyrec, &status))
2852 fits_get_errstatus(status, errmsg);
2853 m_LastError = errmsg;
2854 m_WCSState = Failure;
2858 fits_free_memory(header, &status);
2863 for(
auto &fitsKeyword : m_HeaderRecords)
2866 rec.
append(fitsKeyword.key.leftJustified(8,
' ').toLatin1());
2868 rec.
append(fitsKeyword.value.toByteArray());
2870 rec.
append(fitsKeyword.comment.toLatin1());
2877 if ((status = wcspih(header_str.
data(), nkeyrec, WCSHDR_all, 0, &nreject, &m_nwcs, &m_WCSHandle)) != 0)
2879 wcsvfree(&m_nwcs, &m_WCSHandle);
2880 m_WCSHandle =
nullptr;
2882 m_LastError =
QString(
"wcspih ERROR %1: %2.").
arg(status).
arg(wcshdr_errmsg[status]);
2883 m_WCSState = Failure;
2887 if (m_WCSHandle ==
nullptr)
2889 m_LastError =
i18n(
"No world coordinate systems found.");
2890 m_WCSState = Failure;
2895 if (m_WCSHandle->crpix[0] == 0)
2897 wcsvfree(&m_nwcs, &m_WCSHandle);
2898 m_WCSHandle =
nullptr;
2900 m_LastError =
i18n(
"No world coordinate systems found.");
2901 m_WCSState = Failure;
2906 if ((status = wcsset(m_WCSHandle)) != 0)
2908 wcsvfree(&m_nwcs, &m_WCSHandle);
2909 m_WCSHandle =
nullptr;
2911 m_LastError =
QString(
"wcsset error %1: %2.").
arg(status).
arg(wcs_errmsg[status]);
2912 m_WCSState = Failure;
2916 m_ObjectsSearched =
false;
2917 m_CatObjectsSearched =
false;
2918 m_WCSState = Success;
2921 qCDebug(KSTARS_FITS) <<
"Finished WCS Data processing...";
2931#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
2933 double imgcrd[NWCSFIX], phi, pixcrd[NWCSFIX], theta, worldcrd[NWCSFIX];
2935 if (m_WCSHandle ==
nullptr)
2937 m_LastError =
i18n(
"No world coordinate systems found.");
2944 if ((status = wcss2p(m_WCSHandle, 1, 2, worldcrd, &phi, &theta, imgcrd, pixcrd, stat)) != 0)
2946 m_LastError =
QString(
"wcss2p error %1: %2.").
arg(status).
arg(wcs_errmsg[status]);
2950 wcsImagePoint.
setX(imgcrd[0]);
2951 wcsImagePoint.
setY(imgcrd[1]);
2953 wcsPixelPoint.
setX(pixcrd[0]);
2954 wcsPixelPoint.
setY(pixcrd[1]);
2959 Q_UNUSED(wcsPixelPoint);
2960 Q_UNUSED(wcsImagePoint);
2965bool FITSData::pixelToWCS(
const QPointF &wcsPixelPoint,
SkyPoint &wcsCoord)
2967#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
2969 double imgcrd[NWCSFIX], phi, pixcrd[NWCSFIX], theta, world[NWCSFIX];
2971 if (m_WCSHandle ==
nullptr)
2973 m_LastError =
i18n(
"No world coordinate systems found.");
2977 pixcrd[0] = wcsPixelPoint.
x();
2978 pixcrd[1] = wcsPixelPoint.
y();
2980 if ((status = wcsp2s(m_WCSHandle, 1, 2, pixcrd, imgcrd, &phi, &theta, world, stat)) != 0)
2982 m_LastError =
QString(
"wcsp2s error %1: %2.").
arg(status).
arg(wcs_errmsg[status]);
2987 wcsCoord.
setRA0(world[0] / 15.0);
2993 Q_UNUSED(wcsPixelPoint);
2999#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3000bool FITSData::searchObjects()
3002 return (Options::fitsCatalog() == CAT_SKYMAP) ? searchSkyMapObjects() : searchCatObjects();
3005bool FITSData::searchSkyMapObjects()
3007 if (m_ObjectsSearched)
3010 m_ObjectsSearched =
true;
3015 pixelToWCS(
QPointF(0, 0), startPoint);
3016 pixelToWCS(
QPointF(width() - 1, height() - 1), endPoint);
3018 return findObjectsInImage(startPoint, endPoint);
3021bool FITSData::searchCatObjects()
3023 if (m_CatObjectsSearched)
3026 m_CatObjectsSearched =
true;
3032 if (catROIRadius() > 0)
3038 ok = pixelToWCS(pt, searchCenter);
3040 ok = pixelToWCS(edgePt, searchEdge);
3046 if (getRecordValue(
"DATE-OBS", date))
3049 tsString = tsString.remove(
'\'').trimmed();
3061 num =
new KSNumbers(KStarsData::Instance()->ut().djd());
3074 pt =
QPoint((width() / 2) - 1, (height() / 2) - 1);
3075 ok = pixelToWCS(pt, searchCenter);
3079 double raEdge = searchCenter.
ra0().
Degrees() + (radius / 60.0);
3081 SkyPoint searchEdge(raEdge / 15.0, decEdge);
3083 ok = wcsToPixel(searchEdge, pEdge, edgePoint);
3087 const double radiusPix = std::hypot((pt.
x() - pEdge.
x()), pt.
y() - pEdge.
y());
3088 setCatSearchROI(pt, radiusPix);
3094 qCDebug(KSTARS_FITS) <<
"Unable to process Catalog Object request...";
3097 if (Options::fitsCatalog() == CAT_SIMBAD)
3098 return findSimbadObjectsInImage(searchCenter, radius);
3103bool FITSData::findWCSBounds(
double &minRA,
double &maxRA,
double &minDec,
double &maxDec)
3105 if (m_WCSHandle ==
nullptr)
3107 m_LastError =
i18n(
"No world coordinate systems found.");
3116 auto updateMinMax = [&](
int x,
int y)
3119 double imgcrd[2], phi, pixcrd[2], theta, world[2];
3124 if (wcsp2s(m_WCSHandle, 1, 2, &pixcrd[0], &imgcrd[0], &phi, &theta, &world[0], &stat[0]))
3127 minRA = std::min(minRA, world[0]);
3128 maxRA = std::max(maxRA, world[0]);
3129 minDec = std::min(minDec, world[1]);
3130 maxDec = std::max(maxDec, world[1]);
3134 for (
int y = 0; y < height(); y++)
3137 updateMinMax(width() - 1, y);
3140 for (
int x = 1; x < width() - 1; x++)
3143 updateMinMax(x, height() - 1);
3150 if (wcsToPixel(NCP, pPoint, imagePoint))
3152 if (pPoint.
x() > 0 && pPoint.
x() < width() && pPoint.
y() > 0 && pPoint.
y() < height())
3157 if (wcsToPixel(SCP, pPoint, imagePoint))
3159 if (pPoint.
x() > 0 && pPoint.
x() < width() && pPoint.
y() > 0 && pPoint.
y() < height())
3168#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3171 if (KStarsData::Instance() ==
nullptr)
3179 if (getRecordValue(
"DATE-OBS", date))
3182 tsString = tsString.remove(
'\'').trimmed();
3194 num =
new KSNumbers(KStarsData::Instance()->ut().djd());
3199 m_SkyObjects.
clear();
3204 int type = oneObject->type();
3205 return (type == SkyObject::STAR || type == SkyObject::PLANET || type == SkyObject::ASTEROID ||
3206 type == SkyObject::COMET || type == SkyObject::SUPERNOVA || type == SkyObject::MOON ||
3207 type == SkyObject::SATELLITE);
3210 double world[2], phi, theta, imgcrd[2], pixcrd[2];
3212 for (
auto &
object : list)
3214 world[0] =
object->ra0().Degrees();
3215 world[1] =
object->dec0().Degrees();
3217 if (wcss2p(m_WCSHandle, 1, 2, &world[0], &phi, &theta, &imgcrd[0], &pixcrd[0], &stat[0]) == 0)
3222 if (x > 0 && y > 0 && x < w && y < h)
3223 m_SkyObjects.
append(
new FITSSkyObject(
object, x, y));
3236bool FITSData::findSimbadObjectsInImage(
SkyPoint searchCenter,
double radius)
3238 m_CatObjects.
clear();
3239 m_CatUpdateTable =
false;
3242 QUrl simbadURL =
QUrl(
"https://simbad.cds.unistra.fr/simbad/sim-coo");
3255 simbadQuery.addQueryItem(
"Coord", coord);
3256 simbadQuery.addQueryItem(
"Radius", radiusStr);
3257 simbadQuery.addQueryItem(
"Radius.unit",
"arcsec");
3258 simbadQuery.addQueryItem(
"CooFrame",
"ICRS");
3259 simbadQuery.addQueryItem(
"CooEpoch", epoch);
3260 simbadQuery.addQueryItem(
"output.format",
"ASCII");
3261 simbadQuery.addQueryItem(
"output.max",
QString(
"%1").arg(10000));
3262 simbadQuery.addQueryItem(
"list.otypesel",
"on");
3263 simbadQuery.addQueryItem(
"otypedisp",
"3");
3264 simbadQuery.addQueryItem(
"list.coo1",
"on");
3265 simbadQuery.addQueryItem(
"frame1",
"ICRS");
3266 simbadQuery.addQueryItem(
"epoch1",
"J2000");
3267 simbadQuery.addQueryItem(
"equi1", epoch);
3268 simbadQuery.addQueryItem(
"list.spsel",
"off");
3269 simbadQuery.addQueryItem(
"list.sizesel",
"on");
3270 simbadQuery.addQueryItem(
"list.bibsel",
"off");
3273 qCDebug(KSTARS_FITS) <<
"Simbad query:" << simbadURL;
3280 qCDebug(KSTARS_FITS) <<
"Error:" << response->
errorString() <<
" occured querying SIMBAD with" << simbadURL;
3281 m_CatQueryInProgress =
false;
3282 emit catalogQueryFailed(
i18n(
"Error querying Simbad"));
3286 m_CatQueryInProgress =
true;
3289 m_CatQueryTimer.
start(60 * 1000);
3290 emit loadingCatalogData();
3295#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3296void FITSData::catTimeout()
3298 m_CatQueryInProgress =
false;
3300 qCDebug(KSTARS_FITS) << text;
3301 emit catalogQueryFailed(text);
3381#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3382void FITSData::simbadResponseReady(
QNetworkReply * simbadReply)
3384 m_CatQueryTimer.
stop();
3385 m_CatQueryInProgress =
false;
3388 qCDebug(KSTARS_FITS) <<
"Error:" << simbadReply->
errorString() <<
" occured in SIMBAD query reply";
3389 m_CatQueryInProgress =
false;
3390 emit catalogQueryFailed(
i18n(
"Error querying Simbad"));
3394 auto data = simbadReply->
readAll();
3401 double dist = 0, magnitude = 0;
3405 for (
int i = 0; i < replyLines.
size(); i++)
3408 if (i == 0 && replyLines[i].contains(
"No astronomical object found : "))
3416 if (replyData.
size() == 9 && replyData[0] ==
"coord" && replyData[6] ==
"radius:")
3417 replyStr =
QString(
"%1 %2 %3").
arg(replyData[1]).
arg(replyData[2]).
arg(replyData[7]);
3418 if (replyStr != m_CatObjQuery)
3420 qCDebug(KSTARS_FITS) <<
"Simbad query:" << m_CatObjQuery <<
"Reply:" << replyStr <<
". Ignoring...";
3426 if (replyLines[i].contains(
"Number of objects : "))
3428 else if (replyLines[i].startsWith(
"Object "))
3432 const int firstBreak = replyLines[i].
indexOf(
"---");
3433 const int secondBreak = replyLines[i].
indexOf(
"---", firstBreak + 3);
3434 name = replyLines[i].
mid(7, firstBreak - 7).trimmed();
3435 type = replyLines[i].
mid(firstBreak + 3, secondBreak - firstBreak - 3).trimmed();
3439 qCDebug(KSTARS_FITS) <<
"Bad Simbad Reply, Ignoring...";
3443 else if (numObjs == 1 && i >= 7)
3445 if (replyLines[i].startsWith(
"Coordinates(ICRS"))
3448 if (replyData.
size() >= 8)
3449 coord = replyData[1] +
" " + replyData[2] +
" " + replyData[3] +
" " +
3450 replyData[5] +
" " + replyData[6] +
" " + replyData[7];
3452 else if (replyLines[i].startsWith(
"Flux V :"))
3455 if (replyData.
size() >= 4)
3458 else if (replyLines[i].startsWith(
"Angular size:"))
3459 sizeStr = replyLines[i].
mid(13, replyLines[i].indexOf(
"~") - 13).trimmed();
3460 else if (replyLines[i].startsWith(
"=========================================="))
3463 if (addCatObject(num, name, type, coord, dist, magnitude, sizeStr))
3467 else if (numObjs > 1 && i >= 9)
3472 name =
type = coord = sizeStr =
"";
3473 magnitude = dist = 0.0;
3474 for (
int j = 0; j < replyData.
size(); j++)
3495 if (addCatObject(num, name, type, coord, dist, magnitude, sizeStr))
3499 if (numObjs != ++dataRow)
3500 qCDebug(KSTARS_FITS) <<
"Simbad Header rows:" << numObjs <<
". Data rows:" << dataRow;
3502 m_CatObjectsSearched =
false;
3506 emit loadedCatalogData();
3511bool FITSData::addCatObject(
const int num,
const QString name,
const QString type,
const QString coord,
const double dist,
3512 const double magnitude,
const QString sizeStr)
3518 if (split.
size() != 6)
3520 qCDebug(KSTARS_FITS) <<
"Coordinates for " <<
name <<
"invalid: " << coord;
3526 ok = r.setFromString(raStr,
false);
3528 ok = d.setFromString(decStr,
true);
3531 qCDebug(KSTARS_FITS) <<
"Coordinates for " <<
name <<
"invalid: " << coord;
3539 if (sizeSplit.
size() >= 1)
3544 qCDebug(KSTARS_FITS) <<
"Angular size for " <<
name <<
"invalid: " << sizeStr;
3546 CatObject catObject;
3547 catObject.num = num;
3548 catObject.catType = CAT_SIMBAD;
3549 catObject.name =
name;
3550 catObject.typeCode =
type;
3551 catObject.typeLabel = getCatObjectLabel(type);
3554 catObject.dist = dist;
3555 catObject.magnitude = magnitude;
3556 catObject.size = size;
3557 catObject.highlight =
false;
3558 catObject.show = getCatObjectFilter(type);
3560 double world[2], phi, theta, imgcrd[2], pixcrd[2];
3563 world[0] = r.Degrees();
3564 world[1] = d.Degrees();
3566 int status = wcss2p(m_WCSHandle, 1, 2, &world[0], &phi, &theta, &imgcrd[0], &pixcrd[0], &stat[0]);
3569 qCDebug(KSTARS_FITS) <<
QString(
"wcss2p error processing object %1 %2: %3.").
arg(name).
arg(status).
arg(wcs_errmsg[status]);
3574 double x = pixcrd[0];
3575 double y = pixcrd[1];
3576 if (x > 0 && y > 0 && x < width() && y < height())
3580 m_CatObjects.
append(catObject);
3585void FITSData::setCatObjectsFilters(
const QStringList filters)
3587 m_CatObjectsFilters = filters;
3590#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3591void FITSData::setCatSearchROI(
const QPoint searchCenter,
const int radius)
3593 m_CatROIPt = searchCenter;
3594 m_CatROIRadius = radius;
3595 Q_UNUSED(searchCenter);
3600#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3606 for (
int i = 0; i < MAX_CAT_OBJ_TYPES; i++)
3608 if (code == catObjTypes[i].code)
3610 label = catObjTypes[i].label;
3613 else if (code == catObjTypes[i].candidateCode)
3615 label = catObjTypes[i].label +
"_Cand";
3624#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3625bool FITSData::getCatObjectFilter(
const QString type)
const
3627 if (m_CatObjectsFilters.
isEmpty() || m_CatObjectsFilters.
contains(type))
3633void FITSData::filterCatObjects()
3635#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3636 bool changed =
false;
3637 for (
int i = 0; i < m_CatObjects.
size(); i++)
3639 bool showState = getCatObjectFilter(m_CatObjects[i].typeCode);
3640 if (m_CatObjects[i].show != showState)
3643 m_CatObjects[i].show = showState;
3651bool FITSData::highlightCatObject(
const int hilite,
const int lolite)
3653#if !defined(KSTARS_LITE) && defined(HAVE_WCSLIB)
3654 if (hilite < 0 || hilite >= m_CatObjects.
size())
3658 if (lolite >= 0 && lolite < m_CatObjects.
size())
3659 m_CatObjects[lolite].highlight =
false;
3663 for (
int i = 0; i < m_CatObjects.
size(); i++)
3665 if (m_CatObjects[i].highlight)
3667 m_CatObjects[i].highlight =
false;
3674 m_CatObjects[hilite].highlight =
true;
3679int FITSData::getFlipVCounter()
const
3681 return flipVCounter;
3684void FITSData::setFlipVCounter(
int value)
3686 flipVCounter = value;
3689int FITSData::getFlipHCounter()
const
3691 return flipHCounter;
3694void FITSData::setFlipHCounter(
int value)
3696 flipHCounter = value;
3699int FITSData::getRotCounter()
const
3704void FITSData::setRotCounter(
int value)
3714template <
typename T>
3715bool FITSData::rotFITS(
int rotate,
int mirror)
3719 uint8_t * rotimage =
nullptr;
3724 else if (rotate == 2)
3726 else if (rotate == 3)
3728 else if (rotate < 0)
3729 rotate = rotate + 360;
3731 nx = m_Statistics.width;
3732 ny = m_Statistics.height;
3734 int BBP = m_Statistics.bytesPerPixel;
3737 rotimage =
new uint8_t[m_Statistics.samples_per_channel * m_Statistics.channels * BBP];
3739 if (rotimage ==
nullptr)
3741 qWarning() <<
"Unable to allocate memory for rotated image buffer!";
3745 auto * rotBuffer =
reinterpret_cast<T *
>(rotimage);
3746 auto * buffer =
reinterpret_cast<T *
>(m_ImageBuffer);
3749 if (rotate < 45 && rotate > -45)
3753 for (
int i = 0; i < m_Statistics.channels; i++)
3755 offset = m_Statistics.samples_per_channel * i;
3756 for (x1 = 0; x1 < nx; x1++)
3759 for (y1 = 0; y1 < ny; y1++)
3760 rotBuffer[(y1 * nx) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3764 else if (mirror == 2)
3766 for (
int i = 0; i < m_Statistics.channels; i++)
3768 offset = m_Statistics.samples_per_channel * i;
3769 for (y1 = 0; y1 < ny; y1++)
3772 for (x1 = 0; x1 < nx; x1++)
3773 rotBuffer[(y2 * nx) + x1 + offset] = buffer[(y1 * nx) + x1 + offset];
3779 for (
int i = 0; i < m_Statistics.channels; i++)
3781 offset = m_Statistics.samples_per_channel * i;
3782 for (y1 = 0; y1 < ny; y1++)
3784 for (x1 = 0; x1 < nx; x1++)
3785 rotBuffer[(y1 * nx) + x1 + offset] = buffer[(y1 * nx) + x1 + offset];
3792 else if (rotate >= 45 && rotate < 135)
3796 for (
int i = 0; i < m_Statistics.channels; i++)
3798 offset = m_Statistics.samples_per_channel * i;
3799 for (y1 = 0; y1 < ny; y1++)
3802 for (x1 = 0; x1 < nx; x1++)
3805 rotBuffer[(y2 * ny) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3810 else if (mirror == 2)
3812 for (
int i = 0; i < m_Statistics.channels; i++)
3814 offset = m_Statistics.samples_per_channel * i;
3815 for (y1 = 0; y1 < ny; y1++)
3817 for (x1 = 0; x1 < nx; x1++)
3818 rotBuffer[(x1 * ny) + y1 + offset] = buffer[(y1 * nx) + x1 + offset];
3824 for (
int i = 0; i < m_Statistics.channels; i++)
3826 offset = m_Statistics.samples_per_channel * i;
3827 for (y1 = 0; y1 < ny; y1++)
3830 for (x1 = 0; x1 < nx; x1++)
3833 rotBuffer[(y2 * ny) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3839 m_Statistics.width = ny;
3840 m_Statistics.height = nx;
3844 else if (rotate >= 135 && rotate < 225)
3848 for (
int i = 0; i < m_Statistics.channels; i++)
3850 offset = m_Statistics.samples_per_channel * i;
3851 for (y1 = 0; y1 < ny; y1++)
3854 for (x1 = 0; x1 < nx; x1++)
3855 rotBuffer[(y2 * nx) + x1 + offset] = buffer[(y1 * nx) + x1 + offset];
3859 else if (mirror == 2)
3861 for (
int i = 0; i < m_Statistics.channels; i++)
3863 offset = m_Statistics.samples_per_channel * i;
3864 for (x1 = 0; x1 < nx; x1++)
3867 for (y1 = 0; y1 < ny; y1++)
3868 rotBuffer[(y1 * nx) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3874 for (
int i = 0; i < m_Statistics.channels; i++)
3876 offset = m_Statistics.samples_per_channel * i;
3877 for (y1 = 0; y1 < ny; y1++)
3880 for (x1 = 0; x1 < nx; x1++)
3883 rotBuffer[(y2 * nx) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3891 else if (rotate >= 225 && rotate < 315)
3895 for (
int i = 0; i < m_Statistics.channels; i++)
3897 offset = m_Statistics.samples_per_channel * i;
3898 for (y1 = 0; y1 < ny; y1++)
3900 for (x1 = 0; x1 < nx; x1++)
3901 rotBuffer[(x1 * ny) + y1 + offset] = buffer[(y1 * nx) + x1 + offset];
3905 else if (mirror == 2)
3907 for (
int i = 0; i < m_Statistics.channels; i++)
3909 offset = m_Statistics.samples_per_channel * i;
3910 for (y1 = 0; y1 < ny; y1++)
3913 for (x1 = 0; x1 < nx; x1++)
3916 rotBuffer[(y2 * ny) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3923 for (
int i = 0; i < m_Statistics.channels; i++)
3925 offset = m_Statistics.samples_per_channel * i;
3926 for (y1 = 0; y1 < ny; y1++)
3929 for (x1 = 0; x1 < nx; x1++)
3932 rotBuffer[(y2 * ny) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3938 m_Statistics.width = ny;
3939 m_Statistics.height = nx;
3943 else if (rotate >= 315 && mirror)
3945 for (
int i = 0; i < m_Statistics.channels; i++)
3947 offset = m_Statistics.samples_per_channel * i;
3948 for (y1 = 0; y1 < ny; y1++)
3950 for (x1 = 0; x1 < nx; x1++)
3954 rotBuffer[(y2 * ny) + x2 + offset] = buffer[(y1 * nx) + x1 + offset];
3960 delete[] m_ImageBuffer;
3961 m_ImageBuffer = rotimage;
3966void FITSData::rotWCSFITS(
int angle,
int mirror)
3970 double ctemp1, ctemp2, ctemp3, ctemp4, naxis1, naxis2;
3971 int WCS_DECIMALS = 6;
3973 naxis1 = m_Statistics.width;
3974 naxis2 = m_Statistics.height;
3976 if (fits_read_key_dbl(fptr,
"CD1_1", &ctemp1, comment, &status))
3985 if (!fits_read_key_dbl(fptr,
"CROTA1", &ctemp1, comment, &status))
3986 fits_update_key_dbl(fptr,
"CROTA1", ctemp1 + 90.0, WCS_DECIMALS, comment, &status);
3988 if (!fits_read_key_dbl(fptr,
"CROTA2", &ctemp1, comment, &status))
3989 fits_update_key_dbl(fptr,
"CROTA2", ctemp1 + 90.0, WCS_DECIMALS, comment, &status);
3997 if (!fits_read_key_dbl(fptr,
"CROTA1", &ctemp1, comment, &status))
3998 fits_update_key_dbl(fptr,
"CROTA1", -ctemp1, WCS_DECIMALS, comment, &status);
4000 if (!fits_read_key_dbl(fptr,
"CROTA2", &ctemp1, comment, &status))
4001 fits_update_key_dbl(fptr,
"CROTA2", -ctemp1, WCS_DECIMALS, comment, &status);
4005 if (!fits_read_key_dbl(fptr,
"LTM1_1", &ctemp1, comment, &status))
4006 fits_update_key_dbl(fptr,
"LTM1_1", -ctemp1, WCS_DECIMALS, comment, &status);
4010 if (!fits_read_key_dbl(fptr,
"CD1_1", &ctemp1, comment, &status))
4011 fits_update_key_dbl(fptr,
"CD1_1", -ctemp1, WCS_DECIMALS, comment, &status);
4013 if (!fits_read_key_dbl(fptr,
"CD1_2", &ctemp1, comment, &status))
4014 fits_update_key_dbl(fptr,
"CD1_2", -ctemp1, WCS_DECIMALS, comment, &status);
4016 if (!fits_read_key_dbl(fptr,
"CD2_1", &ctemp1, comment, &status))
4017 fits_update_key_dbl(fptr,
"CD2_1", -ctemp1, WCS_DECIMALS, comment, &status);
4023 if (!fits_read_key_dbl(fptr,
"LTM1_1", &ctemp1, comment, &status))
4027 if (!fits_read_key_dbl(fptr,
"LTM2_2", &ctemp2, comment, &status))
4028 if (ctemp1 == ctemp2)
4033 if (!fits_read_key_dbl(fptr,
"LTV1", <v1, comment, &status))
4034 fits_delete_key(fptr,
"LTV1", &status);
4035 if (!fits_read_key_dbl(fptr,
"LTV2", <v2, comment, &status))
4036 fits_delete_key(fptr,
"LTV2", &status);
4040 if (!fits_read_key_dbl(fptr,
"CRPIX1", &ctemp3, comment, &status))
4041 fits_update_key_dbl(fptr,
"CRPIX1", (ctemp3 - ltv1) / ctemp1, WCS_DECIMALS, comment, &status);
4043 if (!fits_read_key_dbl(fptr,
"CRPIX2", &ctemp3, comment, &status))
4044 fits_update_key_dbl(fptr,
"CRPIX2", (ctemp3 - ltv2) / ctemp1, WCS_DECIMALS, comment, &status);
4048 if (!fits_read_key_dbl(fptr,
"CD1_1", &ctemp3, comment, &status))
4049 fits_update_key_dbl(fptr,
"CD1_1", ctemp3 / ctemp1, WCS_DECIMALS, comment, &status);
4051 if (!fits_read_key_dbl(fptr,
"CD1_2", &ctemp3, comment, &status))
4052 fits_update_key_dbl(fptr,
"CD1_2", ctemp3 / ctemp1, WCS_DECIMALS, comment, &status);
4054 if (!fits_read_key_dbl(fptr,
"CD2_1", &ctemp3, comment, &status))
4055 fits_update_key_dbl(fptr,
"CD2_1", ctemp3 / ctemp1, WCS_DECIMALS, comment, &status);
4057 if (!fits_read_key_dbl(fptr,
"CD2_2", &ctemp3, comment, &status))
4058 fits_update_key_dbl(fptr,
"CD2_2", ctemp3 / ctemp1, WCS_DECIMALS, comment, &status);
4062 fits_delete_key(fptr,
"LTM1_1", &status);
4063 fits_delete_key(fptr,
"LTM1_2", &status);
4071 if (!fits_read_key_dbl(fptr,
"CRPIX1", &ctemp1, comment, &status) &&
4072 !fits_read_key_dbl(fptr,
"CRPIX2", &ctemp2, comment, &status))
4077 fits_update_key_dbl(fptr,
"CRPIX1", naxis1 - ctemp1, WCS_DECIMALS, comment, &status);
4078 else if (angle == 90)
4080 fits_update_key_dbl(fptr,
"CRPIX1", naxis2 - ctemp2, WCS_DECIMALS, comment, &status);
4081 fits_update_key_dbl(fptr,
"CRPIX2", naxis1 - ctemp1, WCS_DECIMALS, comment, &status);
4083 else if (angle == 180)
4085 fits_update_key_dbl(fptr,
"CRPIX1", ctemp1, WCS_DECIMALS, comment, &status);
4086 fits_update_key_dbl(fptr,
"CRPIX2", naxis2 - ctemp2, WCS_DECIMALS, comment, &status);
4088 else if (angle == 270)
4090 fits_update_key_dbl(fptr,
"CRPIX1", ctemp2, WCS_DECIMALS, comment, &status);
4091 fits_update_key_dbl(fptr,
"CRPIX2", ctemp1, WCS_DECIMALS, comment, &status);
4098 fits_update_key_dbl(fptr,
"CRPIX1", naxis2 - ctemp2, WCS_DECIMALS, comment, &status);
4099 fits_update_key_dbl(fptr,
"CRPIX2", ctemp1, WCS_DECIMALS, comment, &status);
4101 else if (angle == 180)
4103 fits_update_key_dbl(fptr,
"CRPIX1", naxis1 - ctemp1, WCS_DECIMALS, comment, &status);
4104 fits_update_key_dbl(fptr,
"CRPIX2", naxis2 - ctemp2, WCS_DECIMALS, comment, &status);
4106 else if (angle == 270)
4108 fits_update_key_dbl(fptr,
"CRPIX1", ctemp2, WCS_DECIMALS, comment, &status);
4109 fits_update_key_dbl(fptr,
"CRPIX2", naxis1 - ctemp1, WCS_DECIMALS, comment, &status);
4117 if (!fits_read_key_dbl(fptr,
"CDELT1", &ctemp1, comment, &status) &&
4118 !fits_read_key_dbl(fptr,
"CDELT2", &ctemp2, comment, &status))
4123 fits_update_key_dbl(fptr,
"CDELT1", -ctemp1, WCS_DECIMALS, comment, &status);
4124 else if (angle == 90)
4126 fits_update_key_dbl(fptr,
"CDELT1", -ctemp2, WCS_DECIMALS, comment, &status);
4127 fits_update_key_dbl(fptr,
"CDELT2", -ctemp1, WCS_DECIMALS, comment, &status);
4129 else if (angle == 180)
4131 fits_update_key_dbl(fptr,
"CDELT1", ctemp1, WCS_DECIMALS, comment, &status);
4132 fits_update_key_dbl(fptr,
"CDELT2", -ctemp2, WCS_DECIMALS, comment, &status);
4134 else if (angle == 270)
4136 fits_update_key_dbl(fptr,
"CDELT1", ctemp2, WCS_DECIMALS, comment, &status);
4137 fits_update_key_dbl(fptr,
"CDELT2", ctemp1, WCS_DECIMALS, comment, &status);
4144 fits_update_key_dbl(fptr,
"CDELT1", -ctemp2, WCS_DECIMALS, comment, &status);
4145 fits_update_key_dbl(fptr,
"CDELT2", ctemp1, WCS_DECIMALS, comment, &status);
4147 else if (angle == 180)
4149 fits_update_key_dbl(fptr,
"CDELT1", -ctemp1, WCS_DECIMALS, comment, &status);
4150 fits_update_key_dbl(fptr,
"CDELT2", -ctemp2, WCS_DECIMALS, comment, &status);
4152 else if (angle == 270)
4154 fits_update_key_dbl(fptr,
"CDELT1", ctemp2, WCS_DECIMALS, comment, &status);
4155 fits_update_key_dbl(fptr,
"CDELT2", -ctemp1, WCS_DECIMALS, comment, &status);
4166 if (!fits_read_key_dbl(fptr,
"CD1_1", &ctemp1, comment, &status))
4168 fits_read_key_dbl(fptr,
"CD1_2", &ctemp2, comment, &status);
4169 fits_read_key_dbl(fptr,
"CD2_1", &ctemp3, comment, &status);
4170 fits_read_key_dbl(fptr,
"CD2_2", &ctemp4, comment, &status);
4176 fits_update_key_dbl(fptr,
"CD1_2", -ctemp2, WCS_DECIMALS, comment, &status);
4177 fits_update_key_dbl(fptr,
"CD2_1", -ctemp3, WCS_DECIMALS, comment, &status);
4179 else if (angle == 90)
4181 fits_update_key_dbl(fptr,
"CD1_1", -ctemp4, WCS_DECIMALS, comment, &status);
4182 fits_update_key_dbl(fptr,
"CD1_2", -ctemp3, WCS_DECIMALS, comment, &status);
4183 fits_update_key_dbl(fptr,
"CD2_1", -ctemp2, WCS_DECIMALS, comment, &status);
4184 fits_update_key_dbl(fptr,
"CD2_2", -ctemp1, WCS_DECIMALS, comment, &status);
4186 else if (angle == 180)
4188 fits_update_key_dbl(fptr,
"CD1_1", ctemp1, WCS_DECIMALS, comment, &status);
4189 fits_update_key_dbl(fptr,
"CD1_2", ctemp2, WCS_DECIMALS, comment, &status);
4190 fits_update_key_dbl(fptr,
"CD2_1", -ctemp3, WCS_DECIMALS, comment, &status);
4191 fits_update_key_dbl(fptr,
"CD2_2", -ctemp4, WCS_DECIMALS, comment, &status);
4193 else if (angle == 270)
4195 fits_update_key_dbl(fptr,
"CD1_1", ctemp4, WCS_DECIMALS, comment, &status);
4196 fits_update_key_dbl(fptr,
"CD1_2", ctemp3, WCS_DECIMALS, comment, &status);
4197 fits_update_key_dbl(fptr,
"CD2_1", ctemp3, WCS_DECIMALS, comment, &status);
4198 fits_update_key_dbl(fptr,
"CD2_2", ctemp1, WCS_DECIMALS, comment, &status);
4205 fits_update_key_dbl(fptr,
"CD1_1", -ctemp4, WCS_DECIMALS, comment, &status);
4206 fits_update_key_dbl(fptr,
"CD1_2", -ctemp3, WCS_DECIMALS, comment, &status);
4207 fits_update_key_dbl(fptr,
"CD2_1", ctemp1, WCS_DECIMALS, comment, &status);
4208 fits_update_key_dbl(fptr,
"CD2_2", ctemp1, WCS_DECIMALS, comment, &status);
4210 else if (angle == 180)
4212 fits_update_key_dbl(fptr,
"CD1_1", -ctemp1, WCS_DECIMALS, comment, &status);
4213 fits_update_key_dbl(fptr,
"CD1_2", -ctemp2, WCS_DECIMALS, comment, &status);
4214 fits_update_key_dbl(fptr,
"CD2_1", -ctemp3, WCS_DECIMALS, comment, &status);
4215 fits_update_key_dbl(fptr,
"CD2_2", -ctemp4, WCS_DECIMALS, comment, &status);
4217 else if (angle == 270)
4219 fits_update_key_dbl(fptr,
"CD1_1", ctemp4, WCS_DECIMALS, comment, &status);
4220 fits_update_key_dbl(fptr,
"CD1_2", ctemp3, WCS_DECIMALS, comment, &status);
4221 fits_update_key_dbl(fptr,
"CD2_1", -ctemp2, WCS_DECIMALS, comment, &status);
4222 fits_update_key_dbl(fptr,
"CD2_2", -ctemp1, WCS_DECIMALS, comment, &status);
4230 if (!fits_read_key_dbl(fptr,
"CO1_1", &ctemp1, comment, &status))
4235 for (i = 1; i < 13; i++)
4237 sprintf(keyword,
"CO1_%d", i);
4238 fits_delete_key(fptr, keyword, &status);
4240 for (i = 1; i < 13; i++)
4242 sprintf(keyword,
"CO2_%d", i);
4243 fits_delete_key(fptr, keyword, &status);
4249uint8_t * FITSData::getWritableImageBuffer()
4251 return m_ImageBuffer;
4254uint8_t
const * FITSData::getImageBuffer()
const
4256 return m_ImageBuffer;
4259void FITSData::setImageBuffer(uint8_t * buffer)
4261 delete[] m_ImageBuffer;
4262 m_ImageBuffer = buffer;
4265bool FITSData::checkDebayer()
4268 char bayerPattern[64], roworder[64];
4271 if (fits_read_keyword(fptr,
"BAYERPAT", bayerPattern,
nullptr, &status))
4274 if (m_Statistics.dataType != TUSHORT && m_Statistics.dataType != TBYTE)
4276 m_LastError =
i18n(
"Only 8 and 16 bits bayered images supported.");
4279 QString pattern(bayerPattern);
4280 pattern = pattern.remove(
'\'').trimmed();
4283 order = order.remove(
'\'').trimmed();
4285 if (order ==
"BOTTOM-UP" && !(m_Statistics.height % 2))
4287 if (pattern ==
"RGGB")
4289 else if (pattern ==
"GBRG")
4291 else if (pattern ==
"GRBG")
4293 else if (pattern ==
"BGGR")
4298 if (pattern ==
"RGGB")
4299 debayerParams.filter = DC1394_COLOR_FILTER_RGGB;
4300 else if (pattern ==
"GBRG")
4301 debayerParams.filter = DC1394_COLOR_FILTER_GBRG;
4302 else if (pattern ==
"GRBG")
4303 debayerParams.filter = DC1394_COLOR_FILTER_GRBG;
4304 else if (pattern ==
"BGGR")
4305 debayerParams.filter = DC1394_COLOR_FILTER_BGGR;
4309 m_LastError =
i18n(
"Unsupported bayer pattern %1.", pattern);
4313 fits_read_key(fptr, TINT,
"XBAYROFF", &debayerParams.offsetX,
nullptr, &status);
4314 fits_read_key(fptr, TINT,
"YBAYROFF", &debayerParams.offsetY,
nullptr, &status);
4316 if (debayerParams.offsetX == 1)
4321 switch (debayerParams.filter)
4323 case DC1394_COLOR_FILTER_RGGB:
4324 debayerParams.filter = DC1394_COLOR_FILTER_GRBG;
4326 case DC1394_COLOR_FILTER_GBRG:
4327 debayerParams.filter = DC1394_COLOR_FILTER_BGGR;
4329 case DC1394_COLOR_FILTER_GRBG:
4330 debayerParams.filter = DC1394_COLOR_FILTER_RGGB;
4332 case DC1394_COLOR_FILTER_BGGR:
4333 debayerParams.filter = DC1394_COLOR_FILTER_GBRG;
4336 debayerParams.offsetX = 0;
4338 if (debayerParams.offsetX != 0 || debayerParams.offsetY > 1 || debayerParams.offsetY < 0)
4340 m_LastError =
i18n(
"Unsupported bayer offsets %1 %2.", debayerParams.offsetX, debayerParams.offsetY);
4349void FITSData::getBayerParams(BayerParams * param)
4351 param->method = debayerParams.method;
4352 param->filter = debayerParams.filter;
4353 param->offsetX = debayerParams.offsetX;
4354 param->offsetY = debayerParams.offsetY;
4357void FITSData::setBayerParams(BayerParams * param)
4359 debayerParams.method = param->method;
4360 debayerParams.filter = param->filter;
4361 debayerParams.offsetX = param->offsetX;
4362 debayerParams.offsetY = param->offsetY;
4365bool FITSData::debayer(
bool reload)
4369 int anynull = 0,
status = 0;
4371 if (fits_read_img(fptr, m_Statistics.dataType, 1, m_Statistics.samples_per_channel,
nullptr, m_ImageBuffer,
4381 switch (m_Statistics.dataType)
4384 return debayer_8bit();
4387 return debayer_16bit();
4394bool FITSData::debayer_8bit()
4396 dc1394error_t error_code;
4398 uint32_t rgb_size = m_Statistics.samples_per_channel * 3 * m_Statistics.bytesPerPixel;
4399 uint8_t * destinationBuffer =
nullptr;
4403 destinationBuffer =
new uint8_t[rgb_size];
4405 catch (
const std::bad_alloc &e)
4407 logOOMError(rgb_size);
4408 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer: %1", e.what());
4412 auto * bayer_source_buffer =
reinterpret_cast<uint8_t *
>(m_ImageBuffer);
4413 auto * bayer_destination_buffer =
reinterpret_cast<uint8_t *
>(destinationBuffer);
4415 if (bayer_destination_buffer ==
nullptr)
4417 logOOMError(rgb_size);
4418 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer.");
4422 int ds1394_height = m_Statistics.height;
4423 auto dc1394_source = bayer_source_buffer;
4425 if (debayerParams.offsetY == 1)
4427 dc1394_source += m_Statistics.width;
4432 error_code = dc1394_bayer_decoding_8bit(dc1394_source, bayer_destination_buffer, m_Statistics.width, ds1394_height,
4433 debayerParams.filter,
4434 debayerParams.method);
4436 if (error_code != DC1394_SUCCESS)
4438 m_LastError =
i18n(
"Debayer failed (%1)", error_code);
4439 m_Statistics.channels = 1;
4440 delete[] destinationBuffer;
4444 if (m_ImageBufferSize != rgb_size)
4446 delete[] m_ImageBuffer;
4449 m_ImageBuffer =
new uint8_t[rgb_size];
4451 catch (
const std::bad_alloc &e)
4453 delete[] destinationBuffer;
4454 logOOMError(rgb_size);
4455 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer: %1", e.what());
4459 m_ImageBufferSize = rgb_size;
4462 auto bayered_buffer =
reinterpret_cast<uint8_t *
>(m_ImageBuffer);
4466 uint8_t * rBuff = bayered_buffer;
4467 uint8_t * gBuff = bayered_buffer + (m_Statistics.width * m_Statistics.height);
4468 uint8_t * bBuff = bayered_buffer + (m_Statistics.width * m_Statistics.height * 2);
4470 int imax = m_Statistics.samples_per_channel * 3 - 3;
4471 for (
int i = 0; i <= imax; i += 3)
4473 *rBuff++ = bayer_destination_buffer[i];
4474 *gBuff++ = bayer_destination_buffer[i + 1];
4475 *bBuff++ = bayer_destination_buffer[i + 2];
4481 m_Statistics.channels = (m_Mode == FITS_NORMAL || m_Mode == FITS_CALIBRATE) ? 3 : 1;
4482 m_Statistics.dataType = TBYTE;
4483 delete[] destinationBuffer;
4487bool FITSData::debayer_16bit()
4489 dc1394error_t error_code;
4491 uint32_t rgb_size = m_Statistics.samples_per_channel * 3 * m_Statistics.bytesPerPixel;
4492 uint8_t *destinationBuffer =
nullptr;
4495 destinationBuffer =
new uint8_t[rgb_size];
4497 catch (
const std::bad_alloc &e)
4499 logOOMError(rgb_size);
4500 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer: %1", e.what());
4504 auto * bayer_source_buffer =
reinterpret_cast<uint16_t *
>(m_ImageBuffer);
4505 auto * bayer_destination_buffer =
reinterpret_cast<uint16_t *
>(destinationBuffer);
4507 if (bayer_destination_buffer ==
nullptr)
4509 logOOMError(rgb_size);
4510 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer.");
4514 int ds1394_height = m_Statistics.height;
4515 auto dc1394_source = bayer_source_buffer;
4517 if (debayerParams.offsetY == 1)
4519 dc1394_source += m_Statistics.width;
4524 error_code = dc1394_bayer_decoding_16bit(dc1394_source, bayer_destination_buffer, m_Statistics.width, ds1394_height,
4525 debayerParams.filter,
4526 debayerParams.method, 16);
4528 if (error_code != DC1394_SUCCESS)
4530 m_LastError =
i18n(
"Debayer failed (%1)");
4531 m_Statistics.channels = 1;
4532 delete[] destinationBuffer;
4536 if (m_ImageBufferSize != rgb_size)
4538 delete[] m_ImageBuffer;
4541 m_ImageBuffer =
new uint8_t[rgb_size];
4543 catch (
const std::bad_alloc &e)
4545 logOOMError(rgb_size);
4546 delete[] destinationBuffer;
4547 m_LastError =
i18n(
"Unable to allocate memory for temporary bayer buffer: %1", e.what());
4551 m_ImageBufferSize = rgb_size;
4554 auto bayered_buffer =
reinterpret_cast<uint16_t *
>(m_ImageBuffer);
4558 uint16_t * rBuff = bayered_buffer;
4559 uint16_t * gBuff = bayered_buffer + (m_Statistics.width * m_Statistics.height);
4560 uint16_t * bBuff = bayered_buffer + (m_Statistics.width * m_Statistics.height * 2);
4562 int imax = m_Statistics.samples_per_channel * 3 - 3;
4563 for (
int i = 0; i <= imax; i += 3)
4565 *rBuff++ = bayer_destination_buffer[i];
4566 *gBuff++ = bayer_destination_buffer[i + 1];
4567 *bBuff++ = bayer_destination_buffer[i + 2];
4570 m_Statistics.channels = (m_Mode == FITS_NORMAL || m_Mode == FITS_CALIBRATE) ? 3 : 1;
4571 m_Statistics.dataType = TUSHORT;
4572 delete[] destinationBuffer;
4576void FITSData::logOOMError(uint32_t requiredMemory)
4578 qCCritical(KSTARS_FITS) <<
"Debayed memory allocation failure. Required Memory:" <<
KFormat().
formatByteSize(requiredMemory)
4579 <<
"Available system memory:" << KSUtils::getAvailableRAM();
4582double FITSData::getADU()
const
4585 for (
int i = 0; i < m_Statistics.channels; i++)
4586 adu += m_Statistics.mean[i];
4588 return (adu /
static_cast<double>(m_Statistics.channels));
4591QString FITSData::getLastError()
const
4596template <
typename T>
4597void FITSData::convertToQImage(
double dataMin,
double dataMax,
double scale,
double zero,
QImage &image)
4599 const auto * buffer =
reinterpret_cast<const T *
>(getImageBuffer());
4600 const T limit = std::numeric_limits<T>::max();
4601 T bMin = dataMin < 0 ? 0 : dataMin;
4602 T bMax = dataMax > limit ? limit : dataMax;
4603 uint16_t w = width();
4604 uint16_t h = height();
4605 uint32_t size = w * h;
4608 if (channels() == 1)
4611 for (
int j = 0; j < h; j++)
4613 unsigned char * scanLine = image.
scanLine(j);
4615 for (
int i = 0; i < w; i++)
4617 val = qBound(bMin, buffer[j * w + i], bMax);
4618 val = val * scale + zero;
4619 scanLine[i] = qBound<unsigned char>(0,
static_cast<uint8_t
>(val), 255);
4625 double rval = 0, gval = 0, bval = 0;
4628 for (
int j = 0; j < h; j++)
4630 auto * scanLine =
reinterpret_cast<QRgb *
>((image.
scanLine(j)));
4632 for (
int i = 0; i < w; i++)
4634 rval = qBound(bMin, buffer[j * w + i], bMax);
4635 gval = qBound(bMin, buffer[j * w + i + size], bMax);
4636 bval = qBound(bMin, buffer[j * w + i + size * 2], bMax);
4638 value = qRgb(rval * scale + zero, gval * scale + zero, bval * scale + zero);
4640 scanLine[i] = value;
4657 if (future.
result() ==
false)
4660 data.getMinMax(&min, &max);
4668 if (data.channels() == 1)
4673 for (
int i = 0; i < 256; i++)
4674 fitsImage.
setColor(i, qRgb(i, i, i));
4681 double dataMin = data.m_Statistics.mean[0] - data.m_Statistics.stddev[0];
4682 double dataMax = data.m_Statistics.mean[0] + data.m_Statistics.stddev[0] * 3;
4684 double bscale = 255. / (dataMax - dataMin);
4685 double bzero = (-dataMin) * (255. / (dataMax - dataMin));
4688 switch (data.m_Statistics.dataType)
4691 data.convertToQImage<uint8_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4695 data.convertToQImage<int16_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4699 data.convertToQImage<uint16_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4703 data.convertToQImage<int32_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4707 data.convertToQImage<uint32_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4711 data.convertToQImage<
float>(dataMin, dataMax, bscale, bzero, fitsImage);
4715 data.convertToQImage<int64_t>(dataMin, dataMax, bscale, bzero, fitsImage);
4719 data.convertToQImage<
double>(dataMin, dataMax, bscale, bzero, fitsImage);
4733 qCCritical(KSTARS_FITS) <<
"Failed to convert" << filename <<
"to FITS since format" << format <<
"is not supported in Qt";
4741 qCCritical(KSTARS_FITS) <<
"Failed to open image" << filename;
4745 output =
QDir(getTemporaryPath()).
filePath(filename +
".fits");
4748 fitsfile *fptr =
nullptr;
4750 long fpixel = 1, naxis = input.
allGray() ? 2 : 3, nelements, exposure;
4751 long naxes[3] = { input.
width(), input.
height(), naxis == 3 ? 3 : 1 };
4752 char error_status[512] = {0};
4754 if (fits_create_file(&fptr,
QString(
'!' + output).toLocal8Bit().data(), &status))
4756 fits_get_errstatus(status, error_status);
4757 qCCritical(KSTARS_FITS) <<
"Failed to create FITS file. Error:" << error_status;
4761 if (fits_create_img(fptr, BYTE_IMG, naxis, naxes, &status))
4763 qCWarning(KSTARS_FITS) <<
"fits_create_img failed:" << error_status;
4765 fits_flush_file(fptr, &status);
4766 fits_close_file(fptr, &status);
4771 fits_update_key(fptr, TLONG,
"EXPOSURE", &exposure,
"Total Exposure Time", &status);
4776 nelements = naxes[0] * naxes[1];
4777 if (fits_write_img(fptr, TBYTE, fpixel, nelements, input.
bits(), &status))
4779 fits_get_errstatus(status, error_status);
4780 qCWarning(KSTARS_FITS) <<
"fits_write_img GRAY failed:" << error_status;
4782 fits_flush_file(fptr, &status);
4783 fits_close_file(fptr, &status);
4790 nelements = naxes[0] * naxes[1] * 3;
4792 uint8_t *srcBuffer = input.
bits();
4794 uint32_t srcBytes = naxes[0] * naxes[1] * 4 - 4;
4796 uint8_t *rgbBuffer =
new uint8_t[nelements];
4797 if (rgbBuffer ==
nullptr)
4799 qCWarning(KSTARS_FITS) <<
"Not enough memory for RGB buffer";
4800 fits_flush_file(fptr, &status);
4801 fits_close_file(fptr, &status);
4805 uint8_t *subR = rgbBuffer;
4806 uint8_t *subG = rgbBuffer + naxes[0] * naxes[1];
4807 uint8_t *subB = rgbBuffer + naxes[0] * naxes[1] * 2;
4808 for (uint32_t i = 0; i < srcBytes; i += 4)
4810 *subB++ = srcBuffer[i];
4811 *subG++ = srcBuffer[i + 1];
4812 *subR++ = srcBuffer[i + 2];
4815 if (fits_write_img(fptr, TBYTE, fpixel, nelements, rgbBuffer, &status))
4817 fits_get_errstatus(status, error_status);
4818 qCWarning(KSTARS_FITS) <<
"fits_write_img RGB failed:" << error_status;
4820 fits_flush_file(fptr, &status);
4821 fits_close_file(fptr, &status);
4822 delete [] rgbBuffer;
4826 delete [] rgbBuffer;
4829 if (fits_flush_file(fptr, &status))
4831 fits_get_errstatus(status, error_status);
4832 qCWarning(KSTARS_FITS) <<
"fits_flush_file failed:" << error_status;
4834 fits_close_file(fptr, &status);
4838 if (fits_close_file(fptr, &status))
4840 fits_get_errstatus(status, error_status);
4841 qCWarning(KSTARS_FITS) <<
"fits_close_file failed:" << error_status;
4848void FITSData::injectWCS(
double orientation,
double ra,
double dec,
double pixscale,
bool eastToTheRight)
4852 updateWCSHeaderData(orientation, ra, dec, pixscale, eastToTheRight);
4856 fits_update_key(fptr, TDOUBLE,
"OBJCTRA", &ra,
"Object RA", &status);
4857 fits_update_key(fptr, TDOUBLE,
"OBJCTDEC", &dec,
"Object DEC", &status);
4861 fits_update_key(fptr, TINT,
"EQUINOX", &epoch,
"Equinox", &status);
4863 fits_update_key(fptr, TDOUBLE,
"CRVAL1", &ra,
"CRVAL1", &status);
4864 fits_update_key(fptr, TDOUBLE,
"CRVAL2", &dec,
"CRVAL1", &status);
4866 char radecsys[8] =
"FK5";
4867 char ctype1[16] =
"RA---TAN";
4868 char ctype2[16] =
"DEC--TAN";
4870 fits_update_key(fptr, TSTRING,
"RADECSYS", radecsys,
"RADECSYS", &status);
4871 fits_update_key(fptr, TSTRING,
"CTYPE1", ctype1,
"CTYPE1", &status);
4872 fits_update_key(fptr, TSTRING,
"CTYPE2", ctype2,
"CTYPE2", &status);
4874 double crpix1 = width() / 2.0;
4875 double crpix2 = height() / 2.0;
4877 fits_update_key(fptr, TDOUBLE,
"CRPIX1", &crpix1,
"CRPIX1", &status);
4878 fits_update_key(fptr, TDOUBLE,
"CRPIX2", &crpix2,
"CRPIX2", &status);
4881 double secpix1 = eastToTheRight ? pixscale : -pixscale;
4882 double secpix2 = pixscale;
4884 fits_update_key(fptr, TDOUBLE,
"SECPIX1", &secpix1,
"SECPIX1", &status);
4885 fits_update_key(fptr, TDOUBLE,
"SECPIX2", &secpix2,
"SECPIX2", &status);
4887 double degpix1 = secpix1 / 3600.0;
4888 double degpix2 = secpix2 / 3600.0;
4890 fits_update_key(fptr, TDOUBLE,
"CDELT1", °pix1,
"CDELT1", &status);
4891 fits_update_key(fptr, TDOUBLE,
"CDELT2", °pix2,
"CDELT2", &status);
4894 double rotation = 360 - orientation;
4898 fits_update_key(fptr, TDOUBLE,
"CROTA1", &rotation,
"CROTA1", &status);
4899 fits_update_key(fptr, TDOUBLE,
"CROTA2", &rotation,
"CROTA2", &status);
4902 emit headerChanged();
4907void FITSData::updateWCSHeaderData(
const double orientation,
const double ra,
const double dec,
const double pixscale,
4908 const bool eastToTheRight)
4910 updateRecordValue(
"OBJCTRA", ra,
"Object RA");
4911 updateRecordValue(
"OBJCTDEC", dec,
"Object DEC");
4912 updateRecordValue(
"EQUINOX", 2000,
"Equinox");
4913 updateRecordValue(
"CRVAL1", ra,
"CRVAL1");
4914 updateRecordValue(
"CRVAL2", dec,
"CRVAL2");
4916 updateRecordValue(
"RADECSYS",
"'FK5'",
"RADECSYS");
4917 updateRecordValue(
"CTYPE1",
"'RA---TAN'",
"CTYPE1");
4918 updateRecordValue(
"CTYPE2",
"'DEC--TAN'",
"CTYPE2");
4920 updateRecordValue(
"CRPIX1", m_Statistics.width / 2.0,
"CRPIX1");
4921 updateRecordValue(
"CRPIX2", m_Statistics.height / 2.0,
"CRPIX2");
4923 double secpix1 = eastToTheRight ? pixscale : -pixscale;
4924 double secpix2 = pixscale;
4926 double degpix1 = secpix1 / 3600.0;
4927 double degpix2 = secpix2 / 3600.0;
4928 updateRecordValue(
"CDELT1", degpix1,
"CDELT1");
4929 updateRecordValue(
"CDELT2", degpix2,
"CDELT2");
4932 double rotation = 360 - orientation;
4935 updateRecordValue(
"CROTA1", rotation,
"CROTA1");
4936 updateRecordValue(
"CROTA2", rotation,
"CROTA2");
4939bool FITSData::contains(
const QPointF &point)
const
4941 return (point.
x() >= 0 && point.
y() >= 0 && point.
x() <= m_Statistics.width && point.
y() <= m_Statistics.height);
4944void FITSData::saveStatistics(FITSImage::Statistic &other)
4946 other = m_Statistics;
4949void FITSData::restoreStatistics(FITSImage::Statistic &other)
4951 m_Statistics = other;
4956void FITSData::constructHistogram()
4958 switch (m_Statistics.dataType)
4961 constructHistogramInternal<uint8_t>();
4965 constructHistogramInternal<int16_t>();
4969 constructHistogramInternal<uint16_t>();
4973 constructHistogramInternal<int32_t>();
4977 constructHistogramInternal<uint32_t>();
4981 constructHistogramInternal<float>();
4985 constructHistogramInternal<int64_t>();
4989 constructHistogramInternal<double>();
4997template <
typename T> int32_t FITSData::histogramBinInternal(T value,
int channel)
const
4999 return qMax(
static_cast<T
>(0), qMin(
static_cast<T
>(m_HistogramBinCount),
5000 static_cast<T
>(rint((value - m_Statistics.min[channel]) / m_HistogramBinWidth[channel]))));
5003template <
typename T> int32_t FITSData::histogramBinInternal(
int x,
int y,
int channel)
const
5005 if (!m_ImageBuffer || !isHistogramConstructed())
5007 uint32_t samples = m_Statistics.width * m_Statistics.height;
5008 uint32_t offset = channel * samples;
5009 auto *
const buffer =
reinterpret_cast<T
const *
>(m_ImageBuffer);
5010 int index = y * m_Statistics.width + x;
5011 const T &sample = buffer[index + offset];
5012 return histogramBinInternal(sample, channel);
5015int32_t FITSData::histogramBin(
int x,
int y,
int channel)
const
5017 switch (m_Statistics.dataType)
5020 return histogramBinInternal<uint8_t>(x, y, channel);
5024 return histogramBinInternal<int16_t>(x, y, channel);
5028 return histogramBinInternal<uint16_t>(x, y, channel);
5032 return histogramBinInternal<int32_t>(x, y, channel);
5036 return histogramBinInternal<uint32_t>(x, y, channel);
5040 return histogramBinInternal<float>(x, y, channel);
5044 return histogramBinInternal<int64_t>(x, y, channel);
5048 return histogramBinInternal<double>(x, y, channel);
5057template <
typename T>
void FITSData::constructHistogramInternal()
5059 auto *
const buffer =
reinterpret_cast<T
const *
>(m_ImageBuffer);
5060 uint32_t samples = m_Statistics.width * m_Statistics.height;
5061 const uint32_t sampleBy = samples > 500000 ? samples / 500000 : 1;
5063 m_HistogramBinCount = qMax(0., qMin(m_Statistics.max[0] - m_Statistics.min[0], 256.0));
5064 if (m_HistogramBinCount <= 0)
5065 m_HistogramBinCount = 256;
5067 for (
int n = 0; n < m_Statistics.channels; n++)
5069 m_HistogramIntensity[n].fill(0, m_HistogramBinCount + 1);
5070 m_HistogramFrequency[n].fill(0, m_HistogramBinCount + 1);
5071 m_CumulativeFrequency[n].fill(0, m_HistogramBinCount + 1);
5073 const double minBinSize = (m_Statistics.max[n] > 1.1) ? 1.0 : .0001;
5074 m_HistogramBinWidth[n] = qMax(minBinSize, (m_Statistics.max[n] - m_Statistics.min[n]) / m_HistogramBinCount);
5079 for (
int n = 0; n < m_Statistics.channels; n++)
5083 for (
int i = 0; i < m_HistogramBinCount; i++)
5084 m_HistogramIntensity[n][i] = m_Statistics.min[n] + (m_HistogramBinWidth[n] * i);
5088 for (
int n = 0; n < m_Statistics.channels; n++)
5092 uint32_t offset = n * samples;
5094 for (uint32_t i = 0; i < samples; i += sampleBy)
5096 int32_t
id = histogramBinInternal<T>(buffer[i + offset], n);
5097 m_HistogramFrequency[n][id] += sampleBy;
5107 for (
int n = 0; n < m_Statistics.channels; n++)
5111 uint32_t accumulator = 0;
5112 for (
int i = 0; i < m_HistogramBinCount; i++)
5114 accumulator += m_HistogramFrequency[n][i];
5115 m_CumulativeFrequency[n].replace(i, accumulator);
5126 if (m_CumulativeFrequency[RED_CHANNEL][m_HistogramBinCount / 4] > 0)
5127 m_JMIndex = m_CumulativeFrequency[RED_CHANNEL][m_HistogramBinCount / 8] /
static_cast<double>
5128 (m_CumulativeFrequency[RED_CHANNEL][m_HistogramBinCount /
5133 qCDebug(KSTARS_FITS) <<
"FITHistogram: JMIndex " << m_JMIndex;
5135 m_HistogramConstructed =
true;
5136 emit histogramReady();
5139void FITSData::recordLastError(
int errorCode)
5141 char fitsErrorMessage[512] = {0};
5142 fits_get_errstatus(errorCode, fitsErrorMessage);
5143 m_LastError = fitsErrorMessage;
5146double FITSData::getAverageMean(
bool roi)
const
5148 const FITSImage::Statistic* ptr = (roi ? &m_ROIStatistics : &m_Statistics);
5149 if (ptr->channels == 1)
5150 return ptr->mean[0];
5152 return (ptr->mean[0] + ptr->mean[1] + ptr->mean[2]) / 3.0;
5155double FITSData::getAverageMedian(
bool roi)
const
5157 const FITSImage::Statistic* ptr = (roi ? &m_ROIStatistics : &m_Statistics);
5158 if (ptr->channels == 1)
5159 return ptr->median[0];
5161 return (ptr->median[0] + ptr->median[1] + ptr->median[2]) / 3.0;
5164double FITSData::getAverageStdDev(
bool roi)
const
5166 const FITSImage::Statistic* ptr = (roi ? &m_ROIStatistics : &m_Statistics);
5167 if (ptr->channels == 1)
5168 return ptr->stddev[0];
5170 return (ptr->stddev[0] + ptr->stddev[1] + ptr->stddev[2]) / 3.0;
There are several time-dependent values used in position calculations, that are not specific to an ob...
SkyMapComposite * skyComposite()
Extension of QDateTime for KStars KStarsDateTime can represent the date/time as a Julian Day,...
double epoch() const
This is (approximately) the year expressed as a floating-point value.
static KStarsDateTime currentDateTimeUtc()
QList< SkyObject * > findObjectsInArea(const SkyPoint &p1, const SkyPoint &p2)
Provides all necessary information about an object in the sky: its coordinates, name(s),...
The sky coordinates of a point in the sky.
const CachingDms & ra0() const
virtual void updateCoordsNow(const KSNumbers *num)
updateCoordsNow Shortcut for updateCoords( const KSNumbers *num, false, nullptr, nullptr,...
dms angularDistanceTo(const SkyPoint *sp, double *const positionAngle=nullptr) const
Computes the angular distance between two SkyObjects.
void setRA0(dms r)
Sets RA0, the catalog Right Ascension.
const CachingDms & dec0() const
void setDec0(dms d)
Sets Dec0, the catalog Declination.
An angle, stored as degrees, but expressible in many ways.
static dms fromString(const QString &s, bool deg)
Static function to create a DMS object from a QString.
const QString toDMSString(const bool forceSign=false, const bool machineReadable=false, const bool highPrecision=false) const
const QString toHMSString(const bool machineReadable=false, const bool highPrecision=false) const
const double & Degrees() const
QString i18n(const char *text, const TYPE &arg...)
Type type(const QSqlDatabase &db)
bool remove(const QString &column, const QVariant &value)
char * toString(const EngineQuery &query)
KIOCORE_EXPORT StatJob * stat(const QUrl &url, JobFlags flags=DefaultFlags)
KIOCORE_EXPORT QString number(KIO::filesize_t size)
KGUIADDONS_EXPORT QWindow * window(QObject *job)
void error(QWidget *parent, const QString &text, const QString &title, const KGuiItem &buttonOk, Options options=Notify)
KIOCORE_EXPORT QStringList list(const QString &fileClass)
QString name(StandardAction id)
QString label(StandardShortcut id)
const QList< QKeySequence > & end()
NETWORKMANAGERQT_EXPORT NetworkManager::Status status()
KOSM_EXPORT double distance(const std::vector< const OSM::Node * > &path, Coordinate coord)
QByteArray & append(QByteArrayView data)
const_iterator constBegin() const const
const char * constData() const const
bool isEmpty() const const
QByteArray leftJustified(qsizetype width, char fill, bool truncate) const const
void push_back(QByteArrayView str)
qsizetype size() const const
QDateTime currentDateTime()
QDateTime fromString(QStringView string, QStringView format, QCalendar cal)
bool isValid() const const
QString filePath(const QString &fileName) const const
QString path() const const
virtual qint64 size() const const override
virtual void close() override
qint64 size() const const
QString suffix() const const
bool isRunning() const const
bool allGray() const const
void fill(Qt::GlobalColor color)
bool load(QIODevice *device, const char *format)
bool loadFromData(QByteArrayView data, const char *format)
bool save(QIODevice *device, const char *format, int quality) const const
void setColor(int index, QRgb colorValue)
void setColorCount(int colorCount)
QString errorString() const const
qint64 write(const QByteArray &data)
void append(QList< T > &&value)
const_reference at(qsizetype i) const const
bool contains(const AT &value) const const
qsizetype count() const const
iterator erase(const_iterator begin, const_iterator end)
bool isEmpty() const const
QList< T > mid(qsizetype pos, qsizetype length) const const
qsizetype size() const const
void finished(QNetworkReply *reply)
NetworkError error() const const
QMetaObject::Connection connect(const QObject *sender, PointerToMemberFunction signal, Functor functor)
void setObjectName(QAnyStringView name)
bool contains(const QPoint &point, bool proper) const const
bool isNull() const const
bool isValid() const const
QPoint topLeft() const const
bool isNull() const const
QString arg(Args &&... args) const const
bool contains(QChar ch, Qt::CaseSensitivity cs) const const
QString fromStdString(const std::string &str)
QString & insert(qsizetype position, QChar ch)
bool isEmpty() const const
QString mid(qsizetype position, qsizetype n) const const
QString & remove(QChar ch, Qt::CaseSensitivity cs)
QString & replace(QChar before, QChar after, Qt::CaseSensitivity cs)
QStringList split(QChar sep, Qt::SplitBehavior behavior, Qt::CaseSensitivity cs) const const
bool startsWith(QChar c, Qt::CaseSensitivity cs) const const
double toDouble(bool *ok) const const
int toInt(bool *ok, int base) const const
QByteArray toLatin1() const const
QByteArray toLocal8Bit() const const
QString trimmed() const const
bool contains(QLatin1StringView str, Qt::CaseSensitivity cs) const const
qsizetype indexOf(const QRegularExpression &re, qsizetype from) const const
QTextStream & center(QTextStream &stream)
QFuture< void > map(Iterator begin, Iterator end, MapFunctor &&function)
QFuture< T > run(Function function,...)
virtual QString fileName() const const override
void setFileTemplate(const QString &name)
void setSingleShot(bool singleShot)
QString query(ComponentFormattingOptions options) const const
void setQuery(const QString &query, ParsingMode mode)
QString toString(StringFormat mode) const const
bool isValid() const const
QDateTime toDateTime() const const
double toDouble(bool *ok) const const
int toInt(bool *ok) const const
QString toString() const const
Object to hold FITS Header records.