13#include <QtConcurrent>
15#include "fitscentroiddetector.h"
16#include "fits_debug.h"
34bool FITSCentroidDetector::checkCollision(Edge * s1, Edge * s2)
const
38 int diff_x = s1->x - s2->x;
39 int diff_y = s1->y - s2->y;
41 dis = std::abs(sqrt(diff_x * diff_x + diff_y * diff_y));
55 FITSImage::Statistic
const &stats = m_ImageData->getStatistics();
56 switch (stats.dataType)
60 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<uint8_t const>, boundary);
63 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<int16_t const>, boundary);
66 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<uint16_t const>, boundary);
69 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<int32_t const>, boundary);
72 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<uint32_t const>, boundary);
75 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<float const>, boundary);
78 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<int64_t const>, boundary);
81 return QtConcurrent::run(
this, &FITSCentroidDetector::findSources<double const>, boundary);
87bool FITSCentroidDetector::findSources(
const QRect &boundary)
89 FITSImage::Statistic
const &stats = m_ImageData->getStatistics();
90 FITSMode
const m_Mode =
static_cast<FITSMode
>(m_ImageData->
property(
"mode").
toInt());
92 int MINIMUM_STDVAR = getValue(
"MINIMUM_STDVAR", 5).
toInt();
93 int minEdgeWidth = getValue(
"MINIMUM_PIXEL_RANGE", 5).
toInt();
94 double JMIndex = getValue(
"JMINDEX", 100.0).
toDouble();
96 int initStdDev = MINIMUM_STDVAR;
97 double threshold = 0, sum = 0, avg = 0, min = 0;
100 int minimumEdgeCount = MINIMUM_EDGE_LIMIT;
102 auto * buffer =
reinterpret_cast<T
const *
>(m_ImageData->getImageBuffer());
104 float dispersion_ratio = 1.5;
108 if (JMIndex < DIFFUSE_THRESHOLD)
110 minEdgeWidth = JMIndex * 35 + 1;
111 minimumEdgeCount = minEdgeWidth - 1;
116 minimumEdgeCount = 4;
119 while (initStdDev >= 1)
124 minEdgeWidth = qMax(3, minEdgeWidth);
125 minimumEdgeCount = qMax(3, minimumEdgeCount);
127 if (JMIndex < DIFFUSE_THRESHOLD)
130 threshold = stats.max[0] - stats.mean[0] * ((MINIMUM_STDVAR - initStdDev) * 0.5 + 1);
133 if (threshold - min < 0)
135 threshold = stats.mean[0] * ((MINIMUM_STDVAR - initStdDev) * 0.5 + 1);
139 dispersion_ratio = 1.4 - (MINIMUM_STDVAR - initStdDev) * 0.08;
143 threshold = stats.mean[0] + stats.stddev[0] * initStdDev * (0.3 - (MINIMUM_STDVAR - initStdDev) * 0.05);
146 dispersion_ratio = 1.8 - (MINIMUM_STDVAR - initStdDev) * 0.2;
149 qCDebug(KSTARS_FITS) <<
"SNR: " << stats.SNR;
150 qCDebug(KSTARS_FITS) <<
"The threshold level is " << threshold <<
"(actual " << threshold - min
151 <<
") minimum edge width" << minEdgeWidth <<
" minimum edge limit " << minimumEdgeCount;
155 int subX, subY, subW, subH;
159 if (m_Mode == FITS_GUIDE || m_Mode == FITS_FOCUS)
162 subX = round(stats.width * 0.15);
163 subY = round(stats.height * 0.15);
164 subW = stats.width - subX;
165 subH = stats.height - subY;
180 subW = subX + boundary.
width();
181 subH = subY + boundary.
height();
185 for (
int i = subY; i < subH; i++)
189 for (
int j = subX; j < subW; j++)
191 pixVal = buffer[j + (i * stats.width)] - min;
194 if (pixVal >= threshold)
204 if (starDiameter >= minEdgeWidth)
206 float center = avg / sum + 0.5;
209 int i_center = std::floor(center);
212 if (((buffer[i_center + (i * stats.width)] - min) /
213 (buffer[i_center + (i * stats.width) - starDiameter / 2] - min) >=
215 ((buffer[i_center + (i * stats.width)] - min) /
216 (buffer[i_center + (i * stats.width) + starDiameter / 2] - min) >=
220 <<
"Edge center is " << buffer[i_center + (i * stats.width)] - min
221 <<
" Edge is " << buffer[i_center + (i * stats.width) - starDiameter / 2] - min
223 << ((buffer[i_center + (i * stats.width)] - min) /
224 (buffer[i_center + (i * stats.width) - starDiameter / 2] - min))
225 <<
" located at X: " << center <<
" Y: " << i + 0.5;
227 auto * newEdge =
new Edge();
230 newEdge->y = i + 0.5;
231 newEdge->scanned = 0;
232 newEdge->val = buffer[i_center + (i * stats.width)] - min;
233 newEdge->width = starDiameter;
243 avg = sum = starDiameter = 0;
248 qCDebug(KSTARS_FITS) <<
"Total number of edges found is: " << edges.
count();
251 if (edges.
count() == 1 && initStdDev > 1)
257 if (edges.
count() >= MAX_EDGE_LIMIT)
259 qCWarning(KSTARS_FITS) <<
"Too many edges, aborting... " << edges.
count();
264 if (edges.
count() >= minimumEdgeCount)
280 auto const greaterThan = [](Edge
const * a, Edge
const * b)
282 return a->sum > b->sum;
284 std::sort(edges.
begin(), edges.
end(), greaterThan);
288 for (
int i = 0; i < edges.
count(); i++)
290 qCDebug(KSTARS_FITS) <<
"# " << i <<
" Edge at (" << edges[i]->x <<
"," << edges[i]->y <<
") With a value of "
291 << edges[i]->val <<
" and width of " << edges[i]->width <<
" pixels. with sum " << edges[i]->sum;
294 if (edges[i]->scanned == 1)
296 qCDebug(KSTARS_FITS) <<
"Skipping check for center " << i <<
" because it was already counted";
300 qCDebug(KSTARS_FITS) <<
"Investigating edge # " << i <<
" now ...";
305 cen_v = edges[i]->sum;
306 cen_w = edges[i]->width;
315 for (
int j = 0; j < edges.
count(); j++)
317 if (edges[j]->scanned)
320 if (checkCollision(edges[j], edges[i]))
322 if (edges[j]->sum >= cen_v)
324 cen_v = edges[j]->sum;
325 cen_w = edges[j]->width;
328 edges[j]->scanned = 1;
331 avg_x += edges[j]->x * edges[j]->val;
332 avg_y += edges[j]->y * edges[j]->val;
333 sum += edges[j]->val;
339 int cen_limit = (MINIMUM_ROWS_PER_CENTER - (MINIMUM_STDVAR - initStdDev));
341 if (edges.
count() < LOW_EDGE_CUTOFF_1)
343 if (edges.
count() < LOW_EDGE_CUTOFF_2)
349 qCDebug(KSTARS_FITS) <<
"center_count: " << cen_count <<
" and initstdDev= " << initStdDev <<
" and limit is "
357 if (cen_count >= cen_limit)
360 auto * rCenter =
new Edge();
362 rCenter->x = avg_x / sum;
363 rCenter->y = avg_y / sum;
364 width_sum += rCenter->width;
365 rCenter->width = cen_w;
367 qCDebug(KSTARS_FITS) <<
"Found a real center with number with (" << rCenter->x <<
"," << rCenter->y <<
")";
374 cen_x = (int)std::floor(rCenter->x);
375 cen_y = (int)std::floor(rCenter->y);
377 if (cen_x < 0 || cen_x > stats.width || cen_y < 0 || cen_y > stats.height)
385 for (
int k = rCenter->width / 2; k >= -(rCenter->width / 2); k--)
387 FSum += buffer[cen_x - k + (cen_y * stats.width)] - min;
395 TF = buffer[cen_y * stats.width + cen_x] - min;
397 int pixelCounter = 1;
400 for (
int k = 1; k < rCenter->width / 2; k++)
404 qCDebug(KSTARS_FITS) <<
"Stopping at TF " << TF <<
" after #" << k <<
" pixels.";
408 TF += buffer[cen_y * stats.width + cen_x + k] - min;
409 TF += buffer[cen_y * stats.width + cen_x - k] - min;
415 rCenter->HFR = pixelCounter * (HF / TF);
419 qCDebug(KSTARS_FITS) <<
"HFR for this center is " << rCenter->HFR <<
" pixels and the total flux is " << FSum;
421 starCenters.
append(rCenter);
425 if (starCenters.
count() > 1 && m_Mode != FITS_FOCUS)
427 float width_avg = (float)width_sum / starCenters.
count();
428 float lsum = 0, sdev = 0;
430 for (
auto ¢er : starCenters)
431 lsum += (
center->width - width_avg) * (
center->width - width_avg);
433 sdev = (std::sqrt(lsum / (starCenters.
count() - 1))) * 4;
436 foreach (Edge * center, starCenters)
444 m_ImageData->setStarCenters(starCenters);
void append(QList< T > &&value)
qsizetype count() const const
bool removeOne(const AT &t)
QVariant property(const char *name) const const
bool isNull() const const
QTextStream & center(QTextStream &stream)
QFuture< T > run(Function function,...)
double toDouble(bool *ok) const const
int toInt(bool *ok) const const