00001
00002
00003
00004
00005
00006
00007 #include "ColorFilter.h"
00008 #include "Correlation.h"
00009 #include "DocumentModelCoords.h"
00010 #include "EngaugeAssert.h"
00011 #include "GridClassifier.h"
00012 #include <iostream>
00013 #include "Logger.h"
00014 #include <QDebug>
00015 #include <QFile>
00016 #include <QImage>
00017 #include "QtToString.h"
00018 #include "Transformation.h"
00019
00020 int GridClassifier::NUM_PIXELS_PER_HISTOGRAM_BINS = 1;
00021 double GridClassifier::PEAK_HALF_WIDTH = 4;
00022 int GridClassifier::MIN_STEP_PIXELS = 4 * GridClassifier::PEAK_HALF_WIDTH;
00023 const QString GNUPLOT_DELIMITER ("\t");
00024
00025
00026
00027 int GridClassifier::BIN_START_UNSHIFTED = GridClassifier::PEAK_HALF_WIDTH;
00028
00029 using namespace std;
00030
00031 GridClassifier::GridClassifier()
00032 {
00033 }
00034
00035 int GridClassifier::binFromCoordinate (double coord,
00036 double coordMin,
00037 double coordMax) const
00038 {
00039 ENGAUGE_ASSERT (coordMin < coordMax);
00040 ENGAUGE_ASSERT (coordMin <= coord);
00041 ENGAUGE_ASSERT (coord <= coordMax);
00042
00043 int bin = 0.5 + (m_numHistogramBins - 1.0) * (coord - coordMin) / (coordMax - coordMin);
00044
00045 return bin;
00046 }
00047
00048 void GridClassifier::classify (bool isGnuplot,
00049 const QPixmap &originalPixmap,
00050 const Transformation &transformation,
00051 int &countX,
00052 double &startX,
00053 double &stepX,
00054 int &countY,
00055 double &startY,
00056 double &stepY)
00057 {
00058 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::classify";
00059
00060 QImage image = originalPixmap.toImage ();
00061
00062 m_numHistogramBins = image.width() / NUM_PIXELS_PER_HISTOGRAM_BINS;
00063
00064 double xMin, xMax, yMin, yMax;
00065 double binStartX, binStepX, binStartY, binStepY;
00066
00067 m_binsX = new double [m_numHistogramBins];
00068 m_binsY = new double [m_numHistogramBins];
00069
00070 computeGraphCoordinateLimits (image,
00071 transformation,
00072 xMin,
00073 xMax,
00074 yMin,
00075 yMax);
00076 initializeHistogramBins ();
00077 populateHistogramBins (image,
00078 transformation,
00079 xMin,
00080 xMax,
00081 yMin,
00082 yMax);
00083 searchStartStepSpace (isGnuplot,
00084 m_binsX,
00085 "x",
00086 xMin,
00087 xMax,
00088 startX,
00089 stepX,
00090 binStartX,
00091 binStepX);
00092 searchStartStepSpace (isGnuplot,
00093 m_binsY,
00094 "y",
00095 yMin,
00096 yMax,
00097 startY,
00098 stepY,
00099 binStartY,
00100 binStepY);
00101 searchCountSpace (m_binsX,
00102 binStartX,
00103 binStepX,
00104 countX);
00105 searchCountSpace (m_binsY,
00106 binStartY,
00107 binStepY,
00108 countY);
00109
00110 delete [] m_binsX;
00111 delete [] m_binsY;
00112 }
00113
00114 void GridClassifier::computeGraphCoordinateLimits (const QImage &image,
00115 const Transformation &transformation,
00116 double &xMin,
00117 double &xMax,
00118 double &yMin,
00119 double &yMax)
00120 {
00121 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::computeGraphCoordinateLimits";
00122
00123
00124
00125
00126 QPointF posGraphTL, posGraphTR, posGraphBL, posGraphBR;
00127 transformation.transformScreenToRawGraph (QPointF (0, 0) , posGraphTL);
00128 transformation.transformScreenToRawGraph (QPointF (image.width(), 0) , posGraphTR);
00129 transformation.transformScreenToRawGraph (QPointF (0, image.height()) , posGraphBL);
00130 transformation.transformScreenToRawGraph (QPointF (image.width(), image.height()), posGraphBR);
00131
00132
00133 if (transformation.modelCoords().coordsType() == COORDS_TYPE_CARTESIAN) {
00134
00135
00136 xMin = qMin (qMin (qMin (posGraphTL.x(), posGraphTR.x()), posGraphBL.x()), posGraphBR.x());
00137 xMax = qMax (qMax (qMax (posGraphTL.x(), posGraphTR.x()), posGraphBL.x()), posGraphBR.x());
00138 yMin = qMin (qMin (qMin (posGraphTL.y(), posGraphTR.y()), posGraphBL.y()), posGraphBR.y());
00139 yMax = qMax (qMax (qMax (posGraphTL.y(), posGraphTR.y()), posGraphBL.y()), posGraphBR.y());
00140
00141 } else {
00142
00143
00144 xMin = 0.0;
00145 xMax = transformation.modelCoords().thetaPeriod();
00146 yMin = transformation.modelCoords().originRadius();
00147 yMax = qMax (qMax (qMax (posGraphTL.y(), posGraphTR.y()), posGraphBL.y()), posGraphBR.y());
00148
00149 }
00150
00151 ENGAUGE_ASSERT (xMin < xMax);
00152 ENGAUGE_ASSERT (yMin < yMax);
00153 }
00154
00155 double GridClassifier::coordinateFromBin (int bin,
00156 double coordMin,
00157 double coordMax) const
00158 {
00159 ENGAUGE_ASSERT (bin < m_numHistogramBins);
00160 ENGAUGE_ASSERT (coordMin < coordMax);
00161
00162 return coordMin + (coordMax - coordMin) * (double) bin / ((double) m_numHistogramBins - 1.0);
00163 }
00164
00165 void GridClassifier::copyVectorToVector (const double from [],
00166 double to []) const
00167 {
00168 for (int bin = 0; bin < m_numHistogramBins; bin++) {
00169 to [bin] = from [bin];
00170 }
00171 }
00172
00173 void GridClassifier::dumpGnuplotCoordinate (const QString &coordinateLabel,
00174 double corr,
00175 const double *bins,
00176 double coordinateMin,
00177 double coordinateMax,
00178 int binStart,
00179 int binStep) const
00180 {
00181 QString filename = QString ("gridclassifier_%1_corr%2_startMax%3_stepMax%4.gnuplot")
00182 .arg (coordinateLabel)
00183 .arg (corr, 8, 'f', 3, '0')
00184 .arg (binStart)
00185 .arg (binStep);
00186
00187 cout << "Writing gnuplot file: " << filename.toLatin1().data() << "\n";
00188
00189 QFile fileDump (filename);
00190 fileDump.open (QIODevice::WriteOnly | QIODevice::Text);
00191 QTextStream strDump (&fileDump);
00192
00193 int bin;
00194
00195
00196 int binCountMax = 0;
00197 for (bin = 0; bin < m_numHistogramBins; bin++) {
00198 if (bins [bin] > binCountMax) {
00199 binCountMax = qMax ((double) binCountMax,
00200 bins [bin]);
00201 }
00202 }
00203
00204
00205 double *picketFence = new double [m_numHistogramBins];
00206 loadPicketFence (picketFence,
00207 binStart,
00208 binStep,
00209 0,
00210 false);
00211
00212
00213 strDump << "bin"
00214 << GNUPLOT_DELIMITER << "coordinate"
00215 << GNUPLOT_DELIMITER << "binCount"
00216 << GNUPLOT_DELIMITER << "startStep"
00217 << GNUPLOT_DELIMITER << "picketFence" << "\n";
00218
00219
00220 for (bin = 0; bin < m_numHistogramBins; bin++) {
00221
00222 double coordinate = coordinateFromBin (bin,
00223 coordinateMin,
00224 coordinateMax);
00225 double startStepValue (((bin - binStart) % binStep == 0) ? 1 : 0);
00226 strDump << bin
00227 << GNUPLOT_DELIMITER << coordinate
00228 << GNUPLOT_DELIMITER << bins [bin]
00229 << GNUPLOT_DELIMITER << binCountMax * startStepValue
00230 << GNUPLOT_DELIMITER << binCountMax * picketFence [bin] << "\n";
00231 }
00232
00233 delete [] picketFence;
00234 }
00235
00236 void GridClassifier::dumpGnuplotCorrelations (const QString &coordinateLabel,
00237 double valueMin,
00238 double valueMax,
00239 const double signalA [],
00240 const double signalB [],
00241 const double correlations [])
00242 {
00243 QString filename = QString ("gridclassifier_%1_correlations.gnuplot")
00244 .arg (coordinateLabel);
00245
00246 cout << "Writing gnuplot file: " << filename.toLatin1().data() << "\n";
00247
00248 QFile fileDump (filename);
00249 fileDump.open (QIODevice::WriteOnly | QIODevice::Text);
00250 QTextStream strDump (&fileDump);
00251
00252 int bin;
00253
00254
00255 double signalAMax = 1, signalBMax = 1, correlationsMax = 1;
00256 for (bin = 0; bin < m_numHistogramBins; bin++) {
00257 if (bin == 0 || signalA [bin] > signalAMax) {
00258 signalAMax = signalA [bin];
00259 }
00260 if (bin == 0 || signalB [bin] > signalBMax) {
00261 signalBMax = signalB [bin];
00262 }
00263 if (bin == 0 || correlations [bin] > correlationsMax) {
00264 correlationsMax = correlations [bin];
00265 }
00266 }
00267
00268
00269 for (int bin = 0; bin < m_numHistogramBins; bin++) {
00270
00271 strDump << coordinateFromBin (bin,
00272 valueMin,
00273 valueMax)
00274 << GNUPLOT_DELIMITER << signalA [bin] / signalAMax
00275 << GNUPLOT_DELIMITER << signalB [bin] / signalBMax
00276 << GNUPLOT_DELIMITER << correlations [bin] / correlationsMax << "\n";
00277 }
00278 }
00279
00280 void GridClassifier::initializeHistogramBins ()
00281 {
00282 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::initializeHistogramBins";
00283
00284 for (int bin = 0; bin < m_numHistogramBins; bin++) {
00285 m_binsX [bin] = 0;
00286 m_binsY [bin] = 0;
00287 }
00288 }
00289
00290 void GridClassifier::loadPicketFence (double picketFence [],
00291 int binStart,
00292 int binStep,
00293 int count,
00294 bool isCount) const
00295 {
00296 const double PEAK_HEIGHT = 1.0;
00297
00298
00299
00300 ENGAUGE_ASSERT (binStart >= PEAK_HALF_WIDTH);
00301 if (!isCount) {
00302 count = 1 + (m_numHistogramBins - binStart - PEAK_HALF_WIDTH) / binStep;
00303 }
00304
00305
00306 int binStartMinusHalfWidth = binStart - PEAK_HALF_WIDTH;
00307 int binStopPlusHalfWidth = (binStart + (count - 1) * binStep) + PEAK_HALF_WIDTH;
00308
00309
00310
00311 double areaUnnormalized = count * PEAK_HEIGHT * PEAK_HALF_WIDTH;
00312 double normalizationOffset = -1.0 * areaUnnormalized / m_numHistogramBins;
00313
00314 for (int bin = 0; bin < m_numHistogramBins; bin++) {
00315
00316
00317
00318 picketFence [bin] = normalizationOffset;
00319
00320 if ((binStartMinusHalfWidth <= bin) &&
00321 (bin <= binStopPlusHalfWidth)) {
00322
00323
00324 int ordinalClosestPeak = (int) ((bin - binStart + binStep / 2) / binStep);
00325 int binClosestPeak = binStart + ordinalClosestPeak * binStep;
00326
00327
00328 int distanceToClosestPeak = qAbs (bin - binClosestPeak);
00329
00330 if (distanceToClosestPeak < PEAK_HALF_WIDTH) {
00331
00332
00333 picketFence [bin] = 1.0 - (double) distanceToClosestPeak / PEAK_HALF_WIDTH + normalizationOffset;
00334
00335 }
00336 }
00337 }
00338 }
00339
00340 void GridClassifier::populateHistogramBins (const QImage &image,
00341 const Transformation &transformation,
00342 double xMin,
00343 double xMax,
00344 double yMin,
00345 double yMax)
00346 {
00347 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::populateHistogramBins";
00348
00349 ColorFilter filter;
00350 QRgb rgbBackground = filter.marginColor (&image);
00351
00352 for (int x = 0; x < image.width(); x++) {
00353 for (int y = 0; y < image.height(); y++) {
00354
00355 QColor pixel = image.pixel (x, y);
00356
00357
00358 if (!filter.colorCompare (rgbBackground,
00359 pixel.rgb ())) {
00360
00361
00362 QPointF posGraph;
00363 transformation.transformScreenToRawGraph (QPointF (x, y), posGraph);
00364
00365 if (transformation.modelCoords().coordsType() == COORDS_TYPE_POLAR) {
00366
00367
00368 while (posGraph.x() < xMin) {
00369 posGraph.setX (posGraph.x() + transformation.modelCoords().thetaPeriod());
00370 }
00371 while (posGraph.x() > xMax) {
00372 posGraph.setX (posGraph.x() - transformation.modelCoords().thetaPeriod());
00373 }
00374 }
00375
00376 int binX = binFromCoordinate (posGraph.x(), xMin, xMax);
00377 int binY = binFromCoordinate (posGraph.y(), yMin, yMax);
00378
00379 ENGAUGE_ASSERT (0 <= binX);
00380 ENGAUGE_ASSERT (0 <= binY);
00381 ENGAUGE_ASSERT (binX < m_numHistogramBins);
00382 ENGAUGE_ASSERT (binY < m_numHistogramBins);
00383
00384
00385 binX = qMin (binX, m_numHistogramBins - 1);
00386 binY = qMin (binY, m_numHistogramBins - 1);
00387
00388 ++m_binsX [binX];
00389 ++m_binsY [binY];
00390 }
00391 }
00392 }
00393 }
00394
00395 void GridClassifier::searchCountSpace (double bins [],
00396 double binStart,
00397 double binStep,
00398 int &countMax)
00399 {
00400 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::searchCountSpace"
00401 << " start=" << binStart
00402 << " step=" << binStep;
00403
00404
00405 Correlation correlation (m_numHistogramBins);
00406 double *picketFence = new double [m_numHistogramBins];
00407 double corr, corrMax;
00408 bool isFirst = true;
00409 int countStop = 1 + (m_numHistogramBins - binStart) / binStep;
00410 for (int count = 2; count <= countStop; count++) {
00411
00412 loadPicketFence (picketFence,
00413 binStart,
00414 binStep,
00415 count,
00416 true);
00417
00418 correlation.correlateWithoutShift (m_numHistogramBins,
00419 bins,
00420 picketFence,
00421 corr);
00422 if (isFirst || (corr > corrMax)) {
00423 countMax = count;
00424 corrMax = corr;
00425 }
00426
00427 isFirst = false;
00428 }
00429
00430 delete [] picketFence;
00431 }
00432
00433 void GridClassifier::searchStartStepSpace (bool isGnuplot,
00434 double bins [],
00435 const QString &coordinateLabel,
00436 double valueMin,
00437 double valueMax,
00438 double &start,
00439 double &step,
00440 double &binStartMax,
00441 double &binStepMax)
00442 {
00443 LOG4CPP_INFO_S ((*mainCat)) << "GridClassifier::searchStartStepSpace";
00444
00445
00446 double *signalA = new double [m_numHistogramBins];
00447 double *signalB = new double [m_numHistogramBins];
00448 double *correlations = new double [m_numHistogramBins];
00449 double *correlationsMax = new double [m_numHistogramBins];
00450
00451
00452 Correlation correlation (m_numHistogramBins);
00453 double *picketFence = new double [m_numHistogramBins];
00454 int binStart;
00455 double corr = 0, corrMax = 0;
00456 bool isFirst = true;
00457
00458
00459
00460
00461
00462
00463 binStartMax = BIN_START_UNSHIFTED + 1;
00464 binStepMax = qMin (MIN_STEP_PIXELS, m_numHistogramBins / 8);
00465 for (int binStep = qMin (MIN_STEP_PIXELS, m_numHistogramBins / 8); binStep < m_numHistogramBins / 4; binStep++) {
00466
00467 loadPicketFence (picketFence,
00468 BIN_START_UNSHIFTED,
00469 binStep,
00470 PEAK_HALF_WIDTH,
00471 false);
00472
00473 correlation.correlateWithShift (m_numHistogramBins,
00474 bins,
00475 picketFence,
00476 binStart,
00477 corr,
00478 correlations);
00479 if (isFirst || (corr > corrMax)) {
00480
00481 int binStartMaxNext = binStart + BIN_START_UNSHIFTED + 1;
00482
00483
00484 if (binStartMaxNext < m_numHistogramBins) {
00485
00486 binStartMax = binStartMaxNext;
00487 binStepMax = binStep;
00488 corrMax = corr;
00489 copyVectorToVector (bins, signalA);
00490 copyVectorToVector (picketFence, signalB);
00491 copyVectorToVector (correlations, correlationsMax);
00492
00493
00494 if (isGnuplot) {
00495
00496 dumpGnuplotCoordinate(coordinateLabel,
00497 corr,
00498 bins,
00499 valueMin,
00500 valueMax,
00501 binStart,
00502 binStep);
00503 }
00504 }
00505 }
00506
00507 isFirst = false;
00508 }
00509
00510
00511 start = coordinateFromBin (binStartMax,
00512 valueMin,
00513 valueMax);
00514 if (binStartMax + binStepMax < m_numHistogramBins) {
00515
00516
00517 double next = coordinateFromBin (binStartMax + binStepMax,
00518 valueMin,
00519 valueMax);
00520 step = next - start;
00521 } else {
00522
00523
00524 double next = coordinateFromBin (m_numHistogramBins - 1,
00525 valueMin,
00526 valueMax);
00527 step = next - start;
00528 }
00529
00530 if (isGnuplot) {
00531 dumpGnuplotCorrelations (coordinateLabel,
00532 valueMin,
00533 valueMax,
00534 signalA,
00535 signalB,
00536 correlationsMax);
00537 }
00538
00539 delete [] signalA;
00540 delete [] signalB;
00541 delete [] correlations;
00542 delete [] correlationsMax;
00543 delete [] picketFence;
00544 }