src/libkstmath/curve.cpp

0001 /***************************************************************************
0002                           vcurve.cpp: holds info for a curve
0003                              -------------------
0004     begin                : Fri Nov 3 2000
0005     copyright            : (C) 2000 by cbn
0006     email                : netterfield@astro.utoronto.ca
0007  ***************************************************************************/
0008
0009 /***************************************************************************
0010  *                                                                         *
0011  *   This program is free software; you can redistribute it and/or modify  *
0012  *   it under the terms of the GNU General Public License as published by  *
0013  *   the Free Software Foundation; either version 2 of the License, or     *
0014  *   (at your option) any later version.                                   *
0015  *                                                                         *
0016  ***************************************************************************/
0017 #include "curve.h"
0018
0019 // includes for KDE
0020
0021 #include <qdebug.h>
0022 #include <QPolygonF>
0023 #include <QXmlStreamWriter>
0024
0025 // application specific includes
0026 #include "dialoglauncher.h"
0027 #include "colorsequence.h"
0028 #include "datacollection.h"
0029 #include "debug.h"
0030 #include "linestyle.h"
0031 #include "math_kst.h"
0032 #include "datavector.h"
0033 #include "ksttimers.h"
0034 #include "objectstore.h"
0035 #include "relationscriptinterface.h"
0036
0037 #include <time.h>
0038 #include <iostream>
0039
0040 // #define DEBUG_VECTOR_CURVE
0041 //#define BENCHMARK
0042
0043 // NOTE: on a modern (eg, sandybridge or later) cpu
0044 // the cpu's branch prediction is so good, this does
0045 // more or less nothing
0046 #ifndef KDE_IS_LIKELY
0047 #if __GNUC__ - 0 >= 3
0048 # define KDE_ISLIKELY( x )    __builtin_expect(!!(x),1)
0049 # define KDE_ISUNLIKELY( x )  __builtin_expect(!!(x),0)
0050 #else
0051 # define KDE_ISLIKELY( x )   ( x )
0052 # define KDE_ISUNLIKELY( x )  ( x )
0053 #endif
0054 #endif
0055
0056 // for painting
0057 #define MAX_NUM_POLYLINES       3
0058
0059 using namespace std;
0060
0061 namespace Kst {
0062
0063 const QString Curve::staticTypeString = "Curve";
0064 const QString Curve::staticTypeTag = "curve";
0065
0066 static const QLatin1String& XVECTOR = QLatin1String("X");
0067 static const QLatin1String& YVECTOR = QLatin1String("Y");
0068 static const QLatin1String& EXVECTOR = QLatin1String("EX");
0069 static const QLatin1String& EYVECTOR = QLatin1String("EY");
0070 static const QLatin1String& EXMINUSVECTOR = QLatin1String("EXMinus");
0071 static const QLatin1String& EYMINUSVECTOR = QLatin1String("EYMinus");
0072
0073 Curve::Curve(ObjectStore *store)
0074 : Relation(store) {
0075   setHasPoints(false);
0076   setHasHead(false);
0077   setHasBars(false);
0078   setHasLines(true);
0079   setLineWidth(1);
0080   setLineStyle(0);
0081   setPointDensity(0);
0082   setPointSize(CURVE_DEFAULT_POINT_SIZE);
0083
0084   MaxX = MinX = MeanX = MaxY = MinY = MeanY = MinPosX = MinPosY = 0;
0085   NS = 0;
0086   _typeString = tr("Curve");
0087   _type = "Curve";
0088   _initializeShortName();
0089   Color = QColor();
0090   HeadColor = QColor();
0091   BarFillColor = QColor();
0092 }
0093
0094 void Curve::_initializeShortName() {
0095   _shortName = 'C'+QString::number(_curvecnum);
0096   if (_curvecnum>max_curvenum)
0097     max_curvenum = _curvecnum;
0098   _curvecnum++;
0099 }
0100
0101
0102 Curve::~Curve() {
0103 }
0104
0105 ScriptInterface* Curve::createScriptInterface() {
0106   return new CurveSI(this);
0107 }
0108
0109
0110 void Curve::internalUpdate() {
0111   Q_ASSERT(myLockStatus() == KstRWLock::WRITELOCKED);
0112
0113   VectorPtr cxV = *_inputVectors.find(XVECTOR);
0114   VectorPtr cyV = *_inputVectors.find(YVECTOR);
0115   if (!cxV || !cyV) {
0116     return;
0117   }
0118
0119   writeLockInputsAndOutputs();
0120
0121   MaxX = cxV->max();
0122   MinX = cxV->min();
0123   MeanX = cxV->mean();
0124   MinPosX = cxV->minPos();
0125
0126   if (MinPosX > MaxX) {
0127     MinPosX = 0;
0128   }
0129   MaxY = cyV->max();
0130   MinY = cyV->min();
0131   MeanY = cyV->mean();
0132   MinPosY = cyV->minPos();
0133
0134   if (MinPosY > MaxY) {
0135     MinPosY = 0;
0136   }
0137
0138   NS = qMax(cxV->length(), cyV->length());
0139
0140   unlockInputsAndOutputs();
0141
0142   _redrawRequired = true;
0143
0144   return;
0145 }
0146
0147
0148 void Curve::point(int i, double &x, double &y) const {
0149   VectorPtr xv = xVector();
0150   if (xv) {
0151     x = xv->interpolate(i, NS);
0152   }
0153   VectorPtr yv = yVector();
0154   if (yv) {
0155     y = yv->interpolate(i, NS);
0156   }
0157 }
0158
0159
0160 void Curve::getEXPoint(int i, double &x, double &y, double &ex) {
0161   VectorPtr xv = xVector();
0162   if (xv) {
0163     x = xv->interpolate(i, NS);
0164   }
0165   VectorPtr yv = yVector();
0166   if (yv) {
0167     y = yv->interpolate(i, NS);
0168   }
0169   VectorPtr exv = xErrorVector();
0170   if (exv) {
0171     ex = exv->interpolate(i, NS);
0172   }
0173 }
0174
0175
0176 void Curve::getEXMinusPoint(int i, double &x, double &y, double &ex) {
0177   VectorPtr xv = xVector();
0178   if (xv) {
0179     x = xv->interpolate(i, NS);
0180   }
0181   VectorPtr yv = yVector();
0182   if (yv) {
0183     y = yv->interpolate(i, NS);
0184   }
0185   VectorPtr exmv = xMinusErrorVector();
0186   if (exmv) {
0187     ex = exmv->interpolate(i, NS);
0188   }
0189 }
0190
0191
0192 void Curve::getEXPoints(int i, double &x, double &y, double &exminus, double &explus) {
0193   VectorPtr xv = xVector();
0194   if (xv) {
0195     x = xv->interpolate(i, NS);
0196   }
0197   VectorPtr yv = yVector();
0198   if (yv) {
0199     y = yv->interpolate(i, NS);
0200   }
0201   VectorPtr exv = xErrorVector();
0202   if (exv) {
0203     explus = exv->interpolate(i, NS);
0204   }
0205   VectorPtr exmv = xMinusErrorVector();
0206   if (exmv) {
0207     exminus = exmv->interpolate(i, NS);
0208   }
0209 }
0210
0211
0212 void Curve::getEYPoint(int i, double &x, double &y, double &ey) {
0213   VectorPtr xv = xVector();
0214   if (xv) {
0215     x = xv->interpolate(i, NS);
0216   }
0217   VectorPtr yv = yVector();
0218   if (yv) {
0219     y = yv->interpolate(i, NS);
0220   }
0221   VectorPtr eyv = yErrorVector();
0222   if (eyv) {
0223     ey = eyv->interpolate(i, NS);
0224   }
0225 }
0226
0227
0228 void Curve::getEYMinusPoint(int i, double &x, double &y, double &ey) {
0229   VectorPtr xv = xVector();
0230   if (xv) {
0231     x = xv->interpolate(i, NS);
0232   }
0233   VectorPtr yv = yVector();
0234   if (yv) {
0235     y = yv->interpolate(i, NS);
0236   }
0237   VectorPtr eyv = yMinusErrorVector();
0238   if (eyv) {
0239     ey = eyv->interpolate(i, NS);
0240   }
0241 }
0242
0243
0244 void Curve::getEYPoints(int i, double &x, double &y, double &eyminus, double &eyplus) {
0245   VectorPtr xv = xVector();
0246   if (xv) {
0247     x = xv->interpolate(i, NS);
0248   }
0249   VectorPtr yv = yVector();
0250   if (yv) {
0251     y = yv->interpolate(i, NS);
0252   }
0253   VectorPtr eyv = yErrorVector();
0254   if (eyv) {
0255     eyplus = eyv->interpolate(i, NS);
0256   }
0257   VectorPtr eymv = yMinusErrorVector();
0258   if (eymv) {
0259     eyminus = eymv->interpolate(i, NS);
0260   }
0261 }
0262
0263
0264 bool Curve::hasXError() const {
0265   return _inputVectors.contains(EXVECTOR);
0266 }
0267
0268
0269 bool Curve::hasYError() const {
0270   return _inputVectors.contains(EYVECTOR);
0271 }
0272
0273
0274 bool Curve::hasXMinusError() const {
0275   return _inputVectors.contains(EXMINUSVECTOR);
0276 }
0277
0278
0279 bool Curve::hasYMinusError() const {
0280   return _inputVectors.contains(EYMINUSVECTOR);
0281 }
0282
0283
0284 void Curve::save(QXmlStreamWriter &s) {
0285   s.writeStartElement(staticTypeTag);
0286   s.writeAttribute("xvector", _inputVectors[XVECTOR]->Name());
0287   s.writeAttribute("yvector", _inputVectors[YVECTOR]->Name());
0288   if (_inputVectors.contains(EXVECTOR)) {
0289     s.writeAttribute("errorxvector", _inputVectors[EXVECTOR]->Name());
0290   }
0291   if (_inputVectors.contains(EYVECTOR)) {
0292     s.writeAttribute("erroryvector", _inputVectors[EYVECTOR]->Name());
0293   }
0294   if (_inputVectors.contains(EXMINUSVECTOR)) {
0295     s.writeAttribute("errorxminusvector", _inputVectors[EXMINUSVECTOR]->Name());
0296   }
0297   if (_inputVectors.contains(EYMINUSVECTOR)) {
0298     s.writeAttribute("erroryminusvector", _inputVectors[EYMINUSVECTOR]->Name());
0299   }
0300   s.writeAttribute("color", Color.name());
0301   s.writeAttribute("alpha", QString::number(Color.alpha()));
0302   s.writeAttribute("headcolor", HeadColor.name());
0303   s.writeAttribute("headalpha", QString::number(HeadColor.alpha()));
0304   s.writeAttribute("barfillcolor", BarFillColor.name());
0305   s.writeAttribute("barfillalpha", QString::number(BarFillColor.alpha()));
0306
0307   s.writeAttribute("haslines", QVariant(HasLines).toString());
0308   s.writeAttribute("linewidth", QString::number(LineWidth));
0309   s.writeAttribute("linestyle", QString::number(LineStyle));
0310
0311   s.writeAttribute("haspoints", QVariant(HasPoints).toString());
0312   s.writeAttribute("pointtype", QString::number(PointType));
0313   s.writeAttribute("pointdensity", QString::number(PointDensity));
0314   s.writeAttribute("pointsize", QString::number(PointSize));
0315
0316   s.writeAttribute("hasbars", QVariant(HasBars).toString());
0317   s.writeAttribute("ignoreautoscale", QVariant(_ignoreAutoScale).toString());
0318
0319   s.writeAttribute("hashead", QVariant(HasHead).toString());
0320   s.writeAttribute("headtype", QString::number(HeadType));
0321
0322   saveNameInfo(s, CURVENUM);
0323
0324   s.writeEndElement();
0325 }
0326
0327
0328 void Curve::setXVector(VectorPtr new_vx) {
0329   if (new_vx) {
0330     _inputVectors[XVECTOR] = new_vx;
0331   } else {
0332     _inputVectors.remove(XVECTOR);
0333   }
0334 }
0335
0336
0337 void Curve::setYVector(VectorPtr new_vy) {
0338   if (new_vy) {
0339     _inputVectors[YVECTOR] = new_vy;
0340   } else {
0341     _inputVectors.remove(YVECTOR);
0342   }
0343 }
0344
0345
0346 void Curve::setXError(VectorPtr new_ex) {
0347   if (new_ex) {
0348     _inputVectors[EXVECTOR] = new_ex;
0349   } else {
0350     _inputVectors.remove(EXVECTOR);
0351   }
0352 }
0353
0354
0355 void Curve::setYError(VectorPtr new_ey) {
0356   if (new_ey) {
0357     _inputVectors[EYVECTOR] = new_ey;
0358   } else {
0359     _inputVectors.remove(EYVECTOR);
0360   }
0361 }
0362
0363
0364 void Curve::setXMinusError(VectorPtr new_ex) {
0365   if (new_ex) {
0366     _inputVectors[EXMINUSVECTOR] = new_ex;
0367   } else {
0368     _inputVectors.remove(EXMINUSVECTOR);
0369   }
0370 }
0371
0372
0373 void Curve::setYMinusError(VectorPtr new_ey) {
0374   if (new_ey) {
0375     _inputVectors[EYMINUSVECTOR] = new_ey;
0376   } else {
0377     _inputVectors.remove(EYMINUSVECTOR);
0378   }
0379 }
0380
0381
0382 // labels for plots
0383 LabelInfo Curve::xLabelInfo() const {
0384   return xVector()->labelInfo();
0385 }
0386
0387
0388 LabelInfo Curve::yLabelInfo() const {
0389   return yVector()->labelInfo();
0390 }
0391
0392
0393 LabelInfo Curve::titleInfo() const {
0394
0395   LabelInfo label_info = yVector()->titleInfo();
0396
0397   if (!_manualDescriptiveName.isEmpty()) {
0398     label_info.name = _manualDescriptiveName;
0399   }
0400
0401   return label_info;
0402 }
0403
0404
0405 QString Curve::propertyString() const {
0406   return tr("%1 vs %2", "a plot of arg1 vs arg2").arg(yVector()->Name()).arg(xVector()->Name());
0407 }
0408
0409
0410 void Curve::showNewDialog() {
0411   DialogLauncher::self()->showCurveDialog();
0412 }
0413
0414
0415 void Curve::showEditDialog() {
0416   DialogLauncher::self()->showCurveDialog(this);
0417 }
0418
0419
0420 int Curve::samplesPerFrame() const {
0421   const DataVector *rvp = qobject_cast<const DataVector*>(_inputVectors[YVECTOR].data());
0422   return rvp ? rvp->samplesPerFrame() : 1;
0423 }
0424
0425
0426 VectorPtr Curve::xVector() const {
0427   return *_inputVectors.find(XVECTOR);
0428 }
0429
0430
0431 VectorPtr Curve::yVector() const {
0432   return *_inputVectors.find(YVECTOR);
0433 }
0434
0435
0436 VectorPtr Curve::xErrorVector() const {
0437   return *_inputVectors.find(EXVECTOR);
0438 }
0439
0440
0441 VectorPtr Curve::yErrorVector() const {
0442   return *_inputVectors.find(EYVECTOR);
0443 }
0444
0445
0446 VectorPtr Curve::xMinusErrorVector() const {
0447   return *_inputVectors.find(EXMINUSVECTOR);
0448 }
0449
0450
0451 VectorPtr Curve::yMinusErrorVector() const {
0452   return *_inputVectors.find(EYMINUSVECTOR);
0453 }
0454
0455
0456 bool Curve::xIsRising() const {
0457   return _inputVectors[XVECTOR]->isRising();
0458 }
0459
0460
0461 inline int indexNearX(double x, VectorPtr& xv, int NS) {
0462   // monotonically rising: we can do a binary search
0463   // should be reasonably fast
0464   if (xv->isRising()) {
0465     int i_top = NS - 1;
0466     int i_bot = 0;
0467
0468     // don't pre-check for x outside of the curve since this is not
0469     // the common case.  It will be correct - just slightly slower...
0470     while (i_bot + 1 < i_top) {
0471       int i0 = (i_top + i_bot)/2;
0472       double rX = xv->interpolate(i0, NS);
0473       if (x < rX) {
0474         i_top = i0;
0475       } else {
0476         i_bot = i0;
0477       }
0478     }
0479     double xt = xv->interpolate(i_top, NS);
0480     double xb = xv->interpolate(i_bot, NS);
0481     if (xt - x < x - xb) {
0482       return i_top;
0483     } else {
0484       return i_bot;
0485     }
0486   } else {
0487     // Oh Oh... not monotonically rising - we have to search the entire curve!
0488     // May be unbearably slow for large vectors
0489     double rX = xv->interpolate(0, NS);
0490     double dx0 = fabs(x - rX);
0491     int i0 = 0;
0492
0493     for (int i = 1; i < NS; ++i) {
0494       rX = xv->interpolate(i, NS);
0495       double dx = fabs(x - rX);
0496       if (dx < dx0) {
0497         dx0 = dx;
0498         i0 = i;
0499       }
0500     }
0501     return i0;
0502   }
0503 }
0504
0505
0506 /** getIndexNearXY: return index of point within (or closest too)
0507     x +- dx which is closest to y **/
0508 int Curve::getIndexNearXY(double x, double dx_per_pix, double y) const {
0509   VectorPtr xv = *_inputVectors.find(XVECTOR);
0510   VectorPtr yv = *_inputVectors.find(YVECTOR);
0511   if (!xv || !yv) {
0512     return 0; // anything better we can do?
0513   }
0514
0515   double xi, yi, dx, dxi, dy, dyi;
0516   bool first = true;
0517   int i,i0, iN, index;
0518   int sc = sampleCount();
0519
0520   if (xv->isRising()) {
0521     iN = i0 = indexNearX(x, xv, NS);
0522
0523     xi = xv->interpolate(i0, NS);
0524     while (i0 > 0 && x-dx_per_pix < xi) {
0525       xi = xv->interpolate(--i0, NS);
0526     }
0527
0528     xi = xv->interpolate(iN, NS);
0529     while (iN < sc-1 && x+dx_per_pix > xi) {
0530       xi = xv->interpolate(++iN, NS);
0531     }
0532   } else {
0533     i0 = 0;
0534     iN = sampleCount()-1;
0535   }
0536
0537   index = i0;
0538   xi = xv->interpolate(index, NS);
0539   yi = yv->interpolate(index, NS);
0540   dx = fabs(x - xi);
0541   dy = fabs(y - yi);
0542
0543   for (i = i0 + 1; i <= iN; ++i) {
0544     xi = xv->interpolate(i, NS);
0545     dxi = fabs(x - xi);
0546     if (dxi < dx_per_pix) {
0547       dx = dxi;
0548       yi = yv->interpolate(i, NS);
0549       dyi = fabs(y - yi);
0550       if (first || dyi < dy) {
0551         first = false;
0552         index = i;
0553         dy = dyi;
0554       }
0555     } else if (dxi < dx) {
0556       dx = dxi;
0557       index = i;
0558     }
0559   }
0560   return index;
0561 }
0562
0563
0564 void Curve::setHasPoints(bool in_HasPoints) {
0565   HasPoints = in_HasPoints;
0566 }
0567
0568
0569 void Curve::setHasHead(bool in_HasHead) {
0570   HasHead = in_HasHead;
0571 }
0572
0573
0574 void Curve::setHasLines(bool in_HasLines) {
0575   HasLines = in_HasLines;
0576 }
0577
0578
0579 void Curve::setHasBars(bool in_HasBars) {
0580   HasBars = in_HasBars;
0581 }
0582
0583
0584 void Curve::setLineWidth(int in_LineWidth) {
0585 #ifdef Q_OS_WIN32
0586   if (in_LineWidth == 0) {
0587     in_LineWidth = 1;
0588   }
0589 #endif
0590   LineWidth = in_LineWidth;
0591 }
0592
0593
0594 void Curve::setLineStyle(int in_LineStyle) {
0595   LineStyle = in_LineStyle;
0596 }
0597
0598
0599 void Curve::setPointDensity(int in_PointDensity) {
0600   PointDensity = in_PointDensity;
0601 }
0602
0603
0604 void Curve::setPointType(int in_PointType) {
0605   PointType = in_PointType;
0606 }
0607
0608
0609 void Curve::setPointSize(double in_PointSize) {
0610   PointSize = in_PointSize;
0611 }
0612
0613
0614 void Curve::setHeadType(int in_HeadType) {
0615   HeadType = in_HeadType;
0616 }
0617
0618
0619 void Curve::setColor(const QColor& new_c) {
0620   Color = new_c;
0621 }
0622
0623 void Curve::setHeadColor(const QColor& new_c) {
0624   HeadColor = new_c;
0625 }
0626
0627 void Curve::setBarFillColor(const QColor &new_c) {
0628   BarFillColor = new_c;
0629 }
0630
0631 double Curve::maxX() const {
0632   if (hasBars() && sampleCount() > 0) {
0633     return MaxX + (MaxX - MinX)/(2*(sampleCount()-1));
0634   }
0635   return MaxX;
0636 }
0637
0638
0639 double Curve::minX() const {
0640   if (hasBars() && sampleCount() > 0) {
0641     return MinX - (MaxX - MinX)/(2*(sampleCount()-1));
0642   }
0643   return MinX;
0644 }
0645
0646
0647 RelationPtr Curve::makeDuplicate() const {
0648   CurvePtr curve = store()->createObject<Curve>();
0649
0650   if (descriptiveNameIsManual()) {
0651     curve->setDescriptiveName(descriptiveName());
0652   }
0653   curve->setXVector(xVector());
0654   curve->setYVector(yVector());
0655   if (hasXError()) {
0656     curve->setXError(xErrorVector());
0657   }
0658   if (hasYError()) {
0659     curve->setYError(yErrorVector());
0660   }
0661   if (hasXMinusError()) {
0662     curve->setXMinusError(xMinusErrorVector());
0663   }
0664   if (hasYMinusError()) {
0665     curve->setYMinusError(yMinusErrorVector());
0666   }
0667
0668   curve->setColor(Color);
0669   curve->setHeadColor(HeadColor);
0670   curve->setBarFillColor(BarFillColor);
0671   curve->setHasPoints(HasPoints);
0672   curve->setHasLines(HasLines);
0673   curve->setHasBars(HasBars);
0674   curve->setHasHead(HasHead);
0675   curve->setLineWidth(LineWidth);
0676   curve->setLineStyle(LineStyle);
0677   curve->setPointType(PointType);
0678   curve->setPointSize(PointSize);
0679   curve->setPointDensity(PointDensity);
0680
0681   curve->writeLock();
0682   curve->registerChange();
0683   curve->unlock();
0684
0685   return RelationPtr(curve);
0686 }
0687
0688 // FIXME: magic numbers (1/4000 and 1/2000) should be based on the reference plot size,
0689 // and not just pulled out of the air by what looks ~good.
0690 // They are currently "about right" for printing to US Letter.
0691 double Curve::pointDim(QRectF w) const {
0692   return qMax(double(1.01), PointSize * ((w.width() + w.height()) * (1.0 / 4000.0)));
0693 }
0694
0695
0696 double Curve::lineDim(const QRectF &R, double linewidth) {
0697   qreal lw = qMax(qreal(1.01), qreal(linewidth)*(R.width()+R.height())*qreal(1.0/2000.0));
0698   // if you ask for a thicker line, make sure it is at least 2 px wide
0699   if (linewidth>1.9) {
0700     lw = qMax(qreal(2.0),lw);
0701   }
0702
0703   return (int(lw+0.5));
0704 }
0705
0706 void tmpPolyLine(const QPolygonF& poly, QPainter *p) {
0707   int n = poly.count();
0708   for (int i=0; i< n-1; i++) {
0709     p->drawLine(poly[i], poly[i+1]);
0710   }
0711 }
0712
0713 void Curve::paintObjects(const CurveRenderContext& context) {
0714 #ifdef BENCHMARK
0715   QTime bench_time;
0716   bench_time.start();
0717 #endif
0718
0719   QPainter *p = context.painter;
0720   double point_dim = pointDim(p->window());
0721   p->save();
0722
0723   p->setRenderHint(QPainter::Antialiasing, context.antialias);
0724   Qt::PenStyle style = Kst::LineStyle[lineStyle()];
0725
0726   if (hasBars()) {
0727
0728     p->setPen(QPen(barFillColor(), _width, style));
0729
0730     foreach(const QRectF& rect, _filledRects) {
0731         p->fillRect(rect, barFillColor());
0732     }
0733   }
0734   p->setPen(QPen(color(), _width, style)); //, Qt::RoundCap, Qt::RoundJoin));
0735
0736   foreach(const QPolygonF& poly, _polygons) {
0737     tmpPolyLine(poly, p);
0738     //p->drawPolyline(poly);
0739   }
0740   //qDebug() << "E:" << bench_time.elapsed() << "ms w: " << p->viewport().width() << "h: " << p->viewport().height() << " a: " << p->viewport().width() * p->viewport().height();
0741   foreach(const QLineF& line, _lines) {
0742     p->drawLine(line);
0743   }
0744   foreach(const QRectF& rect, _rects) {
0745     p->drawRect(rect);
0746   }
0747
0748   foreach(const QPointF& point, _points) {
0749     CurvePointSymbol::draw(PointType, p, point.x(), point.y(), point_dim);
0750   }
0751
0752   if (hasHead() && _head_valid) {
0753     p->setPen(QPen(headColor(), _width, style));
0754     CurvePointSymbol::draw(HeadType, p, _head.x(), _head.y(), point_dim);
0755   }
0756   p->restore();
0757 #ifdef BENCHMARK
0758   std::cout << "curve drawing:" << bench_time.elapsed() << "ms\n";
0759 #endif
0760 }
0761
0762
0763 void Curve::updatePaintObjects(const CurveRenderContext& context) {
0764   _polygons.clear();
0765   _lines.clear();
0766   _points.clear();
0767   _filledRects.clear();
0768   _rects.clear();
0769
0770   VectorPtr xv = *_inputVectors.find(XVECTOR);
0771   VectorPtr yv = *_inputVectors.find(YVECTOR);
0772   if (!xv || !yv) {
0773     return;
0774   }
0775
0776   double Lx = context.Lx, Hx = context.Hx, Ly = context.Ly, Hy = context.Hy;
0777   double m_X = context.m_X, m_Y = context.m_Y;
0778   double b_X = context.b_X, b_Y = context.b_Y;
0779   double XMin = context.XMin, XMax = context.XMax;
0780   bool xLog = context.xLog, yLog = context.yLog;
0781   double xLogBase = context.xLogBase;
0782   double yLogBase = context.yLogBase;
0783   double maxY = 0.0, minY = 0.0;
0784   double rX = 0.0, rY, rEX, rEY;
0785   double X1 = 0.0, Y1 = 0.0;
0786   double X2 = 0.0, Y2 = 0.0;
0787   double last_x1, last_y1;
0788   bool overlap = false;
0789   int i_pt;
0790
0791 #ifdef BENCHMARK
0792   QTime bench_time, benchtmp;
0793   int b_1 = 0, b_2 = 0, b_3 = 0, b_4 = 0;
0794   bench_time.start();
0795   benchtmp.start();
0796   int numberOfLinesDrawn = 0;
0797   int numberOfPointsDrawn = 0;
0798   int numberOfBarsDrawn = 0;
0799 #endif
0800
0801   _width = lineDim(context.painter->window(), lineWidth());
0802
0803   //qDebug() << context.painter->device()->width() << context.painter->device()->logicalDpiX() <<
0804   //            context.painter->device()->width()/context.painter->device()->logicalDpiX();
0805
0806   double errorFlagDim = pointDim(context.painter->window());
0807   if (sampleCount() > 0) {
0808     int i0, iN;
0809
0810     if (xv->isRising()) {
0811       i0 = indexNearX(XMin, xv, NS);
0812       if (i0 > 0) {
0813         --i0;
0814       }
0815       iN = indexNearX(XMax, xv, NS);
0816       if (iN < sampleCount() - 1) {
0817         ++iN;
0818       }
0819     } else {
0820       i0 = 0;
0821       iN = sampleCount() - 1;
0822     }
0823
0824 #ifdef BENCHMARK
0825     clock_t linesStart = clock();
0826 #endif
0827     if (hasLines()) {
0828       QPolygonF points;
0829       points.reserve(MAX_NUM_POLYLINES);
0830
0831       double lastPlottedX = 0;
0832       double lastPlottedY = 0;
0833       int index = 0;
0834       int i0Start = i0;
0835
0836 // optimize - isnan seems expensive, at least in gcc debug mode
0837 //            cachegrind backs this up.
0838 #undef isnan
0839 #define isnan(x) (x != x)
0840       rX = xv->interpolate(i0, NS);
0841       rY = yv->interpolate(i0, NS);
0842       // if invalid point then look backward for the last valid point.
0843       while (i0 > 0 && (isnan(rX) || isnan(rY))) {
0844         --i0;
0845         rX = xv->interpolate(i0, NS);
0846         rY = yv->interpolate(i0, NS);
0847       }
0848
0849       // if invalid point then look forward for the next valid point...
0850       if (isnan(rX) || isnan(rY)) {
0851         i0 = i0Start;
0852         while (i0 < iN && (isnan(rX) || isnan(rY))) {
0853           ++i0;
0854           rX = xv->interpolate(i0, NS);
0855           rY = yv->interpolate(i0, NS);
0856         }
0857       }
0858
0859       if (xLog) {
0860         rX = logXLo(rX, xLogBase);
0861       }
0862       if (yLog) {
0863         rY = logYLo(rY, yLogBase);
0864       }
0865       last_x1 = m_X*rX + b_X;
0866       last_y1 = m_Y*rY + b_Y;
0867
0868       i_pt = i0;
0869
0870       while (i_pt < iN) {
0871         X2 = last_x1;
0872         Y2 = last_y1;
0873
0874         ++i_pt;
0875         rX = xv->interpolate(i_pt, NS);
0876         rY = yv->interpolate(i_pt, NS);
0877         bool foundNan = false;
0878
0879         // if necessary continue looking for the first valid point...
0880         while (i_pt < iN && (isnan(rX) || isnan(rY))) {
0881 #undef isnan
0882           foundNan = true;
0883           ++i_pt;
0884           rX = xv->interpolate(i_pt, NS);
0885           rY = yv->interpolate(i_pt, NS);
0886         }
0887
0888         if (KDE_ISUNLIKELY(foundNan)) {
0889           if (points.size()>0) {
0890             _polygons.append(points);
0891 #ifdef BENCHMARK
0892             ++numberOfLinesDrawn;
0893 #endif
0894           }
0895           points.resize(0);
0896           if (overlap) {
0897             if (X2 >= Lx && X2 <= Hx) {
0898               if (maxY > Hy && minY <= Hy)
0899                 maxY = Hy;
0900               if (minY < Ly && maxY >= Ly)
0901                 minY = Ly;
0902               if (minY >= Ly && minY <= Hy && maxY >= Ly && maxY <= Hy) {
0903 #ifdef BENCHMARK
0904                 ++numberOfLinesDrawn;
0905 #endif
0906                 double fX2 = floor(X2)+0.5;
0907                 _lines.append(QLineF(fX2, minY, fX2, maxY));
0908               }
0909             }
0910             overlap = false;
0911           }
0912         }
0913
0914         if (xLog) {
0915           rX = logXLo(rX, xLogBase);
0916         }
0917         if (yLog) {
0918           rY = logYLo(rY, yLogBase);
0919         }
0920         X1 = m_X*rX + b_X;
0921         Y1 = m_Y*rY + b_Y;
0922         last_x1 = X1;
0923         last_y1 = Y1;
0924
0925         if (KDE_ISLIKELY(!foundNan)) {
0926           if (KDE_ISLIKELY(samePixel(X1, X2))) {
0927             if (KDE_ISLIKELY(overlap)) {
0928               if (KDE_ISUNLIKELY(Y1 > maxY)) {
0929                 maxY = Y1;
0930               }
0931               if (KDE_ISUNLIKELY(Y1 < minY)) {
0932                 minY = Y1;
0933               }
0934             } else {
0935               if (Y1 < Y2) {
0936                 minY = Y1;
0937                 maxY = Y2;
0938               } else {
0939                 maxY = Y1;
0940                 minY = Y2;
0941               }
0942               overlap = true;
0943             }
0944           } else {
0945             if (KDE_ISLIKELY(overlap)) {
0946               if (KDE_ISLIKELY(X2 >= Lx && X2 <= Hx)) {
0947                 if (KDE_ISUNLIKELY(maxY <= Hy && minY >= Ly)) {
0948
0949                   if (points.size()>MAX_NUM_POLYLINES-2) {
0950                     _polygons.append(points);
0951                     points.resize(0);
0952 #ifdef BENCHMARK
0953                     ++numberOfLinesDrawn;
0954 #endif
0955                   }
0956                   double fX2 = floor(X2)+0.5;
0957
0958                   if (KDE_ISUNLIKELY(minY == maxY)) {
0959                     points.append(QPointF(fX2, maxY));
0960                   } else if (KDE_ISUNLIKELY(Y2 == minY)) {
0961                     points.append(QPointF(fX2, maxY));
0962                     points.append(QPointF(fX2, minY));
0963                   } else if (KDE_ISUNLIKELY(Y2 == maxY)) {
0964                     points.append(QPointF(fX2, minY));
0965                     points.append(QPointF(fX2, maxY));
0966                   } else {
0967                     points.append(QPointF(fX2, minY));
0968                     points.append(QPointF(fX2, maxY));
0969                     if (KDE_ISLIKELY(Y2 >= Ly && Y2 <= Hy)) {
0970                       points.append(QPointF(fX2, Y2));
0971                     }
0972                   }
0973                   lastPlottedX = fX2;
0974                   lastPlottedY = Y2;
0975                 } else {
0976                   if (KDE_ISUNLIKELY(maxY > Hy && minY <= Hy)) {
0977                     maxY = Hy;
0978                   }
0979                   if (KDE_ISUNLIKELY(minY < Ly && maxY >= Ly)) {
0980                     minY = Ly;
0981                   }
0982                   if (KDE_ISUNLIKELY(minY >= Ly && minY <= Hy && maxY >= Ly && maxY <= Hy)) {
0983                     if (points.size()>0) {
0984                       _polygons.append(points);
0985                       points.resize(0);
0986 #ifdef BENCHMARK
0987                       ++numberOfLinesDrawn;
0988 #endif
0989                     }
0990 #ifdef BENCHMARK
0991   ++numberOfLinesDrawn;
0992 #endif
0993                     double fX2 = floor(X2)+0.5;
0994                     _lines.append(QLineF(fX2, minY, fX2, maxY));
0995                   }
0996                 }
0997               }
0998
0999               overlap = false;
1000             }
1001
1002             if (KDE_ISLIKELY(!((X1 < Lx && X2 < Lx) || (X1 > Hx && X2 > Hx)))) {
1003
1004               // trim the line to be within the plot...
1005               if (KDE_ISUNLIKELY(isinf(X1))) {
1006                 Y1 = Y2;
1007                 if (X1 > 0.0) {
1008                   X1 = Hx;
1009                 } else {
1010                   X1 = Lx;
1011                 }
1012               }
1013
1014               if (KDE_ISUNLIKELY(isinf(X2))) {
1015                 Y2 = Y1;
1016                 if (X2 > 0.0) {
1017                   X2 = Hx;
1018                 } else {
1019                   X2 = Lx;
1020                 }
1021               }
1022
1023               if (KDE_ISUNLIKELY(isinf(Y1))) {
1024                 X1 = X2;
1025                 if (Y1 > 0.0) {
1026                   Y1 = Hy;
1027                 } else {
1028                   Y1 = Ly;
1029                 }
1030               }
1031
1032               if (KDE_ISUNLIKELY(isinf(Y2))) {
1033                 X2 = X1;
1034                 if (Y2 > 0.0) {
1035                   Y2 = Hy;
1036                 } else {
1037                   Y2 = Ly;
1038                 }
1039               }
1040
1041               if (KDE_ISUNLIKELY(X1 < Lx && X2 > Lx)) {
1042                 Y1 = (Y2 - Y1) / (X2 - X1) * (Lx - X1) + Y1;
1043                 X1 = Lx;
1044               } else if (KDE_ISUNLIKELY(X2 < Lx && X1 > Lx)) {
1045                 Y2 = (Y1 - Y2) / (X1 - X2) * (Lx - X2) + Y2;
1046                 X2 = Lx;
1047               }
1048
1049               if (KDE_ISUNLIKELY(X1 < Hx && X2 > Hx)) {
1050                 Y2 = (Y2 - Y1) / (X2 - X1) * (Hx - X1) + Y1;
1051                 X2 = Hx;
1052               } else if (KDE_ISUNLIKELY(X2 < Hx && X1 > Hx)) {
1053                 Y1 = (Y1 - Y2) / (X1 - X2) * (Hx - X2) + Y2;
1054                 X1 = Hx;
1055               }
1056
1057               if (KDE_ISUNLIKELY(Y1 < Ly && Y2 > Ly)) {
1058                 X1 = (X2 - X1) / (Y2 - Y1) * (Ly - Y1) + X1;
1059                 Y1 = Ly;
1060               } else if (KDE_ISUNLIKELY(Y2 < Ly && Y1 > Ly)) {
1061                 X2 = (X1 - X2) / (Y1 - Y2) * (Ly - Y2) + X2;
1062                 Y2 = Ly;
1063               }
1064
1065               if (KDE_ISUNLIKELY(Y1 < Hy && Y2 > Hy)) {
1066                 X2 = (X2 - X1) / (Y2 - Y1) * (Hy - Y1) + X1;
1067                 Y2 = Hy;
1068               } else if (KDE_ISUNLIKELY(Y2 < Hy && Y1 > Hy)) {
1069                 X1 = (X1 - X2) / (Y1 - Y2) * (Hy - Y2) + X2;
1070                 Y1 = Hy;
1071               }
1072
1073
1074 #ifdef DEBUG_VECTOR_CURVE
1075 bool xInBounds = (X1 >= Lx && X1 <= Hx && X2 >= Lx && X2 <= Hx);
1076 bool yInBounds = (Y1 >= Ly && Y1 <= Hy && Y2 >= Ly && Y2 <= Hy);
1077
1078 if (!xInBounds)
1079 qDebug() << "x not in bounds"
1080          << "X1:" << X1
1081          << "X2:" << X2
1082          << "Lx:" << Lx
1083          << "Hx:" << Hx
1084          << endl;
1085 if (!yInBounds)
1086 qDebug() << "y not in bounds"
1087          << "Y1:" << Y1
1088          << "Y2:" << Y2
1089          << "Ly:" << Ly
1090          << "Hy:" << Hy
1091          << endl;
1092 #endif
1093
1094               if (X1 >= Lx && X1 <= Hx && X2 >= Lx && X2 <= Hx &&
1095                   Y1 >= Ly && Y1 <= Hy && Y2 >= Ly && Y2 <= Hy) {
1096
1097
1098                 if (KDE_ISUNLIKELY(points.size()==0)) {
1099                     points.append(QPointF(X2, Y2));
1100                     points.append(QPointF(X1, Y1));
1101                 } else if (samePixel(lastPlottedX, X2) &&
1102                     samePixel(lastPlottedY,Y2) &&
1103                     index < MAX_NUM_POLYLINES) {
1104                   points.append(QPointF(X1, Y1));
1105                 } else {
1106                   if (KDE_ISLIKELY(points.size()>1)) {
1107                     _polygons.append(points);
1108                     points.resize(0);
1109 #ifdef BENCHMARK
1110                     ++numberOfLinesDrawn;
1111 #endif
1112                   }
1113                   points.append(QPointF(X2, Y2));
1114                   points.append(QPointF(X1, Y1));
1115                 }
1116                 lastPlottedX = X1;
1117                 lastPlottedY = Y1;
1118               }
1119             }
1120           } // end if (X1 == X2)
1121         } // end if (!foundNan)
1122       } // end while
1123
1124       // we might a have polyline left undrawn...
1125       if (points.size()>1) {
1126         _polygons.append(points);
1127         points.resize(0);
1128 #ifdef BENCHMARK
1129         ++numberOfLinesDrawn;
1130 #endif
1131       }
1132
1133       // we might have some overlapping points still unplotted...
1134       if (overlap) {
1135         if (X2 >= Lx && X2 <= Hx) {
1136           if (maxY > Hy && minY <= Hy) {
1137             maxY = Hy;
1138           }
1139           if (minY < Ly && maxY >= Ly) {
1140             minY = Ly;
1141           }
1142           if (minY >= Ly && minY <= Hy && maxY >= Ly && maxY <= Hy) {
1143 #ifdef BENCHMARK
1144             ++numberOfLinesDrawn;
1145 #endif
1146            _lines.append(QLineF(X2, minY, X2, maxY));
1147           }
1148         }
1149         overlap = false;
1150       }
1151     } // end if hasLines()
1152
1153 #ifdef BENCHMARK
1154     clock_t linesEnd = clock();
1155     std::cout << "\n        Lines clocks: " << (linesEnd - linesStart) << "\n";
1156     b_1 = benchtmp.elapsed();
1157 #endif
1158
1159     VectorPtr exv = _inputVectors.contains(EXVECTOR) ? *_inputVectors.find(EXVECTOR) : 0;
1160     VectorPtr eyv = _inputVectors.contains(EYVECTOR) ? *_inputVectors.find(EYVECTOR) : 0;
1161     VectorPtr exmv = _inputVectors.contains(EXMINUSVECTOR) ? *_inputVectors.find(EXMINUSVECTOR) : 0;
1162     VectorPtr eymv = _inputVectors.contains(EYMINUSVECTOR) ? *_inputVectors.find(EYMINUSVECTOR) : 0;
1163     // draw the bargraph bars, if any...
1164     if (hasBars()) {
1165       bool visible = true;
1166       double rX2 = 0.0;
1167       double drX = 0.0;
1168       QRectF lastRect;
1169
1170       if (!exv) {
1171         // determine the bar position width. NOTE: This is done
1172         //  only if xv->isRising() as in this case the calculation
1173         //  is simple...
1174         drX = (maxX() - minX())/double(sampleCount());
1175         if (xv->isRising()) {
1176           double oldX = 0.0;
1177
1178           for (i_pt = i0; i_pt <= iN; ++i_pt) {
1179             rX = xv->interpolate(i_pt, NS);
1180             if (i_pt > i0) {
1181               if (rX - oldX < drX) {
1182                 drX = rX - oldX;
1183               }
1184             }
1185             oldX = rX;
1186           }
1187         }
1188       }
1189
1190       for (i_pt = i0; i_pt <= iN; ++i_pt) {
1191         visible = true;
1192
1193         if (exv) {
1194           drX = exv->interpolate(i_pt, NS);
1195         }
1196         rX = xv->interpolate(i_pt, NS);
1197         rY = yv->interpolate(i_pt, NS);
1198         rX -= drX/2.0;
1199         rX2 = rX + drX;
1200         if (xLog) {
1201           rX = logXLo(rX, xLogBase);
1202           rX2 = logXLo(rX2, xLogBase);
1203         }
1204         if (yLog) {
1205           rY = logYLo(rY, yLogBase);
1206           Y2 = Hy;
1207         } else {
1208           Y2 = b_Y;
1209         }
1210
1211         X1 = m_X * rX + b_X;
1212         X2 = m_X * rX2 + b_X;
1213         Y1 = m_Y * rY + b_Y;
1214
1215         if ((X1 > Hx) || (X2 < Lx) || (Y1 == Ly && Y2 == Ly) || (Y1 == Hy && Y2 == Hy)) {
1216           visible = false;
1217         }
1218
1219         if (visible) {
1220           QRectF rect(X1, Y1, X2 - X1, Y2 - Y1);
1221           if (!lastRect.contains(rect)) {
1222             //if (barStyle() == 1) { // filled
1223               _filledRects.append(rect);
1224             //} else {
1225               _rects.append(rect);
1226             //}
1227             lastRect = rect;
1228 #ifdef BENCHMARK
1229             ++numberOfBarsDrawn;
1230 #endif
1231           }
1232         }
1233       }
1234     }
1235
1236 #ifdef BENCHMARK
1237     b_2 = benchtmp.elapsed();
1238 #endif
1239
1240     // draw the points, if any...
1241     if (hasPoints()) {
1242       const double w = Hx - Lx;
1243       const double h = Hy - Ly;
1244       int size = 0;
1245       if (hasLines() && pointDensity() != 0) {
1246         // high density: 91 points
1247         // med density 27
1248         // low density 9 points
1249         size = w / int(pow(3.0, POINTDENSITY_MAXTYPE - pointDensity()+1));
1250       }
1251
1252       QRectF rect(Lx, Ly, w, h);
1253       QPointF pt, lastPt;
1254
1255       for (i_pt = i0; i_pt <= iN; ++i_pt) {
1256         rX = xv->interpolate(i_pt, NS);
1257         rY = yv->interpolate(i_pt, NS);
1258         if (xLog) {
1259           rX = logXLo(rX, xLogBase);
1260         }
1261         if (yLog) {
1262           rY = logYLo(rY, yLogBase);
1263         }
1264
1265         pt.setX(m_X * rX + b_X);
1266         pt.setY(m_Y * rY + b_Y);
1267         if (rect.contains(pt) && pt != lastPt &&
1268             (lastPt.isNull() || (fabs(pt.x() - lastPt.x()) > size) || ((size==0) && (fabs(pt.y() - lastPt.y()) > 0)))) {
1269 #ifdef BENCHMARK
1270           ++numberOfPointsDrawn;
1271 #endif
1272             lastPt = pt;
1273             _points.append(pt);
1274         }
1275       }
1276     }
1277
1278     _head_valid = false;
1279     if (hasHead()) {
1280       const double w = Hx - Lx;
1281       const double h = Hy - Ly;
1282
1283       QRectF rect(Lx, Ly, w, h);
1284
1285       rX = xv->interpolate(NS-1, NS);
1286       rY = yv->interpolate(NS-1, NS);
1287       if (xLog) {
1288         rX = logXLo(rX, xLogBase);
1289       }
1290       if (yLog) {
1291         rY = logYLo(rY, yLogBase);
1292       }
1293
1294       _head.setX(m_X * rX + b_X);
1295       _head.setY(m_Y * rY + b_Y);
1296       if (rect.contains(_head)) {
1297         _head_valid = true;
1298       }
1299     }
1300
1301
1302 #ifdef BENCHMARK
1303     b_3 = benchtmp.elapsed();
1304 #endif
1305
1306     // draw the x-errors, if any...
1307     if ((exv || exmv) && !hasBars()) {
1308       double rX1;
1309       double rX2;
1310       bool do_low_flag = true;
1311       bool do_high_flag = true;
1312       bool errorSame = false;
1313
1314       if (exv && exmv) {
1315         if (exv->shortName() == exmv->shortName()) {
1316           errorSame = true;
1317         }
1318       }
1319
1320       for (i_pt = i0; i_pt <= iN; ++i_pt) {
1321         do_low_flag = true;
1322         do_high_flag = true;
1323
1324         rX = xv->interpolate(i_pt, NS);
1325         rY = yv->interpolate(i_pt, NS);
1326         if (errorSame) {
1327           rEX = fabs(exv->interpolate(i_pt, NS));
1328           if (xLog) {
1329             rX1 = logXLo(rX - rEX, xLogBase);
1330             rX2 = logXLo(rX + rEX, xLogBase);
1331           } else {
1332             rX1 = rX - rEX;
1333             rX2 = rX + rEX;
1334           }
1335         } else if (exv && exmv) {
1336           double rEXHi = fabs(exv->interpolate(i_pt, NS));
1337           double rEXLo = fabs(exmv->interpolate(i_pt, NS));
1338           if (xLog) {
1339             rX1 = logXLo(rX - rEXLo, xLogBase);
1340             rX2 = logXLo(rX + rEXHi, xLogBase);
1341           } else {
1342             rX1 = rX - rEXLo;
1343             rX2 = rX + rEXHi;
1344           }
1345         } else if (exv) {
1346           rEX = exv->interpolate(i_pt, NS);
1347           if (xLog) {
1348             rX1 = logXLo(rX, xLogBase);
1349             rX2 = logXLo(rX + fabs(rEX), xLogBase);
1350           } else {
1351             rX1 = rX;
1352             rX2 = rX + fabs(rEX);
1353           }
1354           do_low_flag = false;
1355         } else {
1356           rEX = fabs(exmv->interpolate(i_pt, NS));
1357           if (xLog) {
1358             rX1 = logXLo(rX - rEX, xLogBase);
1359             rX2 = logXLo(rX, xLogBase);
1360           } else {
1361             rX1 = rX - rEX;
1362             rX2 = rX;
1363           }
1364           do_high_flag = false;
1365         }
1366
1367         if (yLog) {
1368           rY = logYLo(rY, yLogBase);
1369         }
1370
1371         X1 = m_X * rX1 + b_X;
1372         X2 = m_X * rX2 + b_X;
1373         Y1 = m_Y * rY + b_Y;
1374
1375         if (X1 < Lx && X2 > Lx) {
1376           X1 = Lx;
1377           do_low_flag = false;
1378         }
1379         if (X1 < Hx && X2 > Hx) {
1380           X2 = Hx;
1381           do_high_flag = false;
1382         }
1383
1384         if (X1 >= Lx && X2 <= Hx && Y1 >= Ly && Y1 <= Hy) {
1385           _lines.append(QLineF(X1, Y1, X2, Y1));
1386           if (do_low_flag) {
1387             _lines.append(QLineF(X1, Y1 + errorFlagDim, X1, Y1 - errorFlagDim));
1388           }
1389           if (do_high_flag) {
1390             _lines.append(QLineF(X2, Y1 + errorFlagDim, X2, Y1 - errorFlagDim));
1391           }
1392         }
1393       }
1394     }
1395
1396     // draw the y-errors, if any...
1397     if (eyv || eymv) {
1398       double rY1;
1399       double rY2;
1400       bool do_low_flag = true;
1401       bool do_high_flag = true;
1402       bool errorSame = false;
1403
1404       if (eyv && eymv) {
1405         if (eyv->shortName() == eymv->shortName()) {
1406           errorSame = true;
1407         }
1408       }
1409
1410       for (i_pt = i0; i_pt <= iN; ++i_pt) {
1411         do_low_flag = true;
1412         do_high_flag = true;
1413
1414         rX = xv->interpolate(i_pt, NS);
1415         rY = yv->interpolate(i_pt, NS);
1416         if (errorSame) {
1417           rEY = eyv->interpolate(i_pt, NS);
1418           if (yLog) {
1419             rY1 = logYLo(rY-fabs(rEY), yLogBase);
1420             rY2 = logYLo(rY+fabs(rEY), yLogBase);
1421           } else {
1422             rY1 = rY-fabs(rEY);
1423             rY2 = rY+fabs(rEY);
1424           }
1425         } else if (eyv && eymv) {
1426           double rEYHi = fabs(eyv->interpolate(i_pt, NS));
1427           double rEYLo = fabs(eymv->interpolate(i_pt, NS));
1428           if (yLog) {
1429             rY1 = logYLo(rY - rEYLo, yLogBase);
1430             rY2 = logYLo(rY + rEYHi, yLogBase);
1431           } else {
1432             rY1 = rY - rEYLo;
1433             rY2 = rY + rEYHi;
1434           }
1435         } else if (eyv) {
1436           rEY = fabs(eyv->interpolate(i_pt, NS));
1437           if (yLog) {
1438             rY1 = logYLo(rY, yLogBase);
1439             rY2 = logYLo(rY + rEY, yLogBase);
1440           } else {
1441             rY1 = rY;
1442             rY2 = rY + rEY;
1443           }
1444           do_low_flag = false;
1445         } else {
1446           rEY = fabs(eymv->interpolate(i_pt, NS));
1447           if (yLog) {
1448             rY1 = logYLo(rY - rEY, yLogBase);
1449             rY2 = logYLo(rY, yLogBase);
1450           } else {
1451             rY1 = rY - rEY;
1452             rY2 = rY;
1453           }
1454           do_high_flag = false;
1455         }
1456
1457         if (xLog) {
1458           rX = logXLo(rX, xLogBase);
1459         }
1460
1461         X1 = m_X * rX + b_X;
1462         Y1 = m_Y * rY1 + b_Y;
1463         Y2 = m_Y * rY2 + b_Y;
1464
1465         if (Y1 < Ly && Y2 > Ly) {
1466           Y1 = Ly;
1467           do_low_flag = false;
1468         }
1469         if (Y1 < Hy && Y2 > Hy) {
1470           Y2 = Hy;
1471           do_high_flag = false;
1472         }
1473
1474         if (X1 >= Lx && X1 <= Hx && Y1 >= Ly && Y2 <= Hy) {
1475           _lines.append(QLineF(X1, Y1, X1, Y2));
1476           if (do_low_flag) {
1477           _lines.append(QLineF(X1 + errorFlagDim, Y1, X1 - errorFlagDim, Y1));
1478           }
1479           if (do_high_flag) {
1480           _lines.append(QLineF(X1 + errorFlagDim, Y2, X1 - errorFlagDim, Y2));
1481           }
1482         }
1483       }
1484     } // end if (hasYError())
1485   } // end if (sampleCount() > 0)
1486
1487 #ifdef BENCHMARK
1488   b_4 = benchtmp.elapsed();
1489 #endif
1490
1491 #ifdef BENCHMARK
1492   int i = bench_time.elapsed();
1493   //qDebug() << endl << "Plotting curve " << (void *)this << ": " << i << "ms";
1494   //qDebug() << "    Without locks: " << b_4 << "ms";
1495   //qDebug() << "    Number of lines drawn:" << numberOfLinesDrawn;
1496   //qDebug() << "    Number of points drawn:" << numberOfPointsDrawn;
1497   //qDebug() << "    Number of bars drawn:" << numberOfBarsDrawn;
1498   std::cout << "\nPlotting curve " << (void *)this << ": " << i << " ms";
1499   std::cout << "\n    Without locks: " << b_4 << "ms";
1500   std::cout << "\n    Number of lines drawn:" << numberOfLinesDrawn;
1501   std::cout << "\n    Number of points drawn:" << numberOfPointsDrawn;
1502   std::cout << "\n    Number of bars drawn:" << numberOfBarsDrawn;
1503   if (b_1 > 0)       std::cout << "\n            Lines: " << b_1 << " ms\n";
1504   if (b_2 - b_1 > 0) std::cout << "\n             Bars: " << (b_2 - b_1) << " ms\n";
1505   if (b_3 - b_2 > 0) std::cout << "\n           Points: " << (b_3 - b_2) << " ms\n";
1506   if (b_4 - b_3 > 0) std::cout << "\n           Errors: " << (b_4 - b_3) << " ms\n";
1507 #endif
1508 }
1509
1510
1511 void Curve::yRange(double xFrom, double xTo, double* yMin, double* yMax) {
1512   if (!yMin || !yMax) {
1513     return;
1514   }
1515
1516   VectorPtr xv = *_inputVectors.find(XVECTOR);
1517   VectorPtr yv = *_inputVectors.find(YVECTOR);
1518   if (!xv || !yv) {
1519     *yMin = *yMax = 0;
1520     return;
1521   }
1522
1523   // get range of the curve to search for min/max
1524   int i0, iN;
1525   if (xv->isRising()) {
1526     i0 = indexNearX(xFrom, xv, NS);
1527     iN = indexNearX(xTo, xv, NS);
1528   } else {
1529     i0 = 0;
1530     iN = sampleCount() - 1;
1531   }
1532   // search for min/max
1533   bool first = true;
1534   double newYMax = 0, newYMin = 0;
1535   for (int i_pt = i0; i_pt <= iN; ++i_pt) {
1536     double rX = xv->interpolate(i_pt, NS);
1537     double rY = yv->interpolate(i_pt, NS);
1538     // make sure this point is visible
1539     if (rX >= xFrom && rX <= xTo) {
1540       // update min/max
1541       if (first || rY > newYMax) {
1542         newYMax = rY;
1543       }
1544       if (first || rY < newYMin) {
1545         newYMin = rY;
1546       }
1547       first = false;
1548     }
1549   }
1550   *yMin = newYMin;
1551   *yMax = newYMax;
1552 }
1553
1554
1555 double Curve::distanceToPoint(double xpos, double dx, double ypos) const {
1556 // find the y distance between the curve and a point. return 1.0E300 if this distance is undefined. i don't want to use -1 because it will make the code which uses this function messy.
1557   VectorPtr xv = *_inputVectors.find(XVECTOR);
1558   if (!xv) {
1559     return 1.0E300; // anything better we can do?
1560   }
1561
1562   double distance = 1.0E300;
1563
1564   int i_near_x = getIndexNearXY(xpos, dx, ypos);
1565   double near_x, near_y;
1566   point(i_near_x, near_x, near_y);
1567
1568   if (fabs(near_x - xpos) < dx) {
1569     distance = fabs(ypos - near_y); // initial estimate.
1570   }
1571
1572   if (hasLines() && xv->isRising()) {
1573     // if hasLines then we should find the distance between the curve and the point, not the data and
1574     //  the point. if isRising because it is (probably) to slow to use this technique if the data is
1575     //  unordered. borrowed from indexNearX. use binary search to find the indices immediately above
1576     //  and below our xpos.
1577     int i_top = NS - 1;
1578     int i_bot = 0;
1579
1580     while (i_bot + 1 < i_top) {
1581       int i0 = (i_top + i_bot)/2;
1582
1583       double rX = xv->interpolate(i0, NS);
1584       if (xpos < rX) {
1585         i_top = i0;
1586       } else {
1587         i_bot = i0;
1588       }
1589     }
1590     // end borrowed
1591
1592     double x_bot, y_bot, x_top, y_top;
1593     point(i_bot, x_bot, y_bot);
1594     point(i_top, x_top, y_top);
1595
1596     if (x_bot <= xpos && x_top >= xpos) {
1597       near_y = (y_top - y_bot) / (x_top - x_bot) * (xpos - x_bot) + y_bot; // calculate y value for line segment between x_bot and x_top.
1598
1599       if (fabs(ypos - near_y) < distance) {
1600         distance = fabs(ypos - near_y);
1601       }
1602     }
1603   }
1604
1605   return distance;
1606 }
1607
1608
1609 QSize Curve::legendSymbolSize(QPainter *p) {
1610   return QSize(p->fontMetrics().height()*3.5, p->fontMetrics().height());
1611 }
1612
1613
1614 void Curve::paintLegendSymbol(QPainter *p, const QSize &size) {
1615   QRect bound(QPoint(0,0),size);
1616
1617   int width;
1618
1619   if (lineWidth() == 0) {
1620     width = 1;
1621   } else {
1622     width = lineWidth();
1623   }
1624
1625   p->save();
1626   if (hasLines()) {
1627     // draw a line from left to right centered vertically
1628     p->setPen(QPen(color(), width, Kst::LineStyle[lineStyle()]));
1629     p->drawLine(QLineF(bound.left(), bound.top() + bound.height()*.5,
1630                 bound.right(), bound.top() + bound.height()*.5));
1631   }
1632   if (hasPoints()) {
1633     // draw a point in the middle
1634     p->setPen(QPen(color(), width));
1635     double point_dim = pointDim(p->window());
1636     CurvePointSymbol::draw(PointType, p, bound.left() + bound.width()*.5, bound.top() + bound.height()*.5, point_dim);
1637   }
1638   p->restore();
1639 }
1640
1641 QString Curve::_automaticDescriptiveName() const {
1642   return tr("%1 vs %2").arg(yVector()->descriptiveName()).arg(xVector()->descriptiveName());
1643 }
1644
1645 QString Curve::descriptionTip() const {
1646   QString tip;
1647
1648   tip = tr("Curve: %1\nX: %2\nY: %3").arg(Name()).arg(xVector()->descriptionTip()).arg(yVector()->descriptionTip());
1649
1650   if (hasXError()) {
1651     tip += tr("\nX+ Error: %1").arg(xErrorVector()->Name());
1652   }
1653
1654   if (hasXMinusError()) {
1655     tip += tr("\nX- Error: %1").arg(xMinusErrorVector()->Name());
1656   }
1657
1658   if (hasYError()) {
1659     tip += tr("\nY+ Error: %1").arg(yErrorVector()->Name());
1660   }
1661
1662   if (hasYMinusError()) {
1663     tip += tr("\nY- Error: %1").arg(yMinusErrorVector()->Name());
1664   }
1665
1666   if (hasLines()) {
1667     tip += tr("\nLines: Width %1 and Style %2").arg(lineWidth()).arg(lineStyle());
1668   }
1669
1670   if (hasPoints()) {
1671     tip += tr("\nPoints: Style %1").arg(pointType());
1672   }
1673
1674   if (hasBars()) {
1675     tip += tr("\nBars");
1676   }
1677
1678   return tip;
1679 }
1680
1681
1682 double Curve::ns_maxX(int zoom_level) const {
1683   return xVector()->ns_max(zoom_level);
1684 }
1685
1686 double Curve::ns_minX(int zoom_level) const {
1687   return xVector()->ns_min(zoom_level);
1688 }
1689
1690 double Curve::ns_maxY(int zoom_level) const {
1691   return yVector()->ns_max(zoom_level);
1692 }
1693
1694 double Curve::ns_minY(int zoom_level) const {
1695   return yVector()->ns_min(zoom_level);
1696 }
1697
1698
1699 }
1700 // vim: ts=2 sw=2 et