labplot  2.8.2
About: LabPlot is an application for plotting and analysis of 2D and 3D functions and data. It is a complete rewrite of LabPlot1 and lacks in the first release a lot of features available in the predecessor. On the other hand, the GUI and the usability is more superior.
  Fossies Dox: labplot-2.8.2.tar.gz  ("unofficial" and yet experimental doxygen-generated source code documentation)  

XYInterpolationCurve.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  File : XYInterpolationCurve.cpp
3  Project : LabPlot
4  Description : A xy-curve defined by an interpolation
5  --------------------------------------------------------------------
6  Copyright : (C) 2016 Stefan Gerlach (stefan.gerlach@uni.kn)
7  Copyright : (C) 2016-2017 Alexander Semke (alexander.semke@web.de)
8 
9  ***************************************************************************/
10 
11 /***************************************************************************
12  * *
13  * This program is free software; you can redistribute it and/or modify *
14  * it under the terms of the GNU General Public License as published by *
15  * the Free Software Foundation; either version 2 of the License, or *
16  * (at your option) any later version. *
17  * *
18  * This program is distributed in the hope that it will be useful, *
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
21  * GNU General Public License for more details. *
22  * *
23  * You should have received a copy of the GNU General Public License *
24  * along with this program; if not, write to the Free Software *
25  * Foundation, Inc., 51 Franklin Street, Fifth Floor, *
26  * Boston, MA 02110-1301 USA *
27  * *
28  ***************************************************************************/
29 
30 /*!
31  \class XYInterpolationCurve
32  \brief A xy-curve defined by an interpolation
33 
34  \ingroup worksheet
35 */
36 
37 #include "XYInterpolationCurve.h"
42 #include "backend/lib/macros.h"
43 #include "backend/gsl/errors.h"
44 
45 extern "C" {
46 #include <gsl/gsl_interp.h>
47 #include <gsl/gsl_spline.h>
48 #include "backend/nsl/nsl_diff.h"
49 #include "backend/nsl/nsl_int.h"
50 }
51 
52 #include <QElapsedTimer>
53 #include <QThreadPool>
54 #include <QIcon>
55 
58 }
59 
62 }
63 
64 //no need to delete the d-pointer here - it inherits from QGraphicsItem
65 //and is deleted during the cleanup in QGraphicsScene
67 
70  d->recalculate();
71 }
72 
73 /*!
74  Returns an icon to be used in the project explorer.
75 */
77  return QIcon::fromTheme("labplot-xy-interpolation-curve");
78 }
79 
80 //##############################################################################
81 //########################## getter methods ##################################
82 //##############################################################################
84 
85 const XYInterpolationCurve::InterpolationResult& XYInterpolationCurve::interpolationResult() const {
86  Q_D(const XYInterpolationCurve);
87  return d->interpolationResult;
88 }
89 
90 //##############################################################################
91 //################# setter methods and undo commands ##########################
92 //##############################################################################
93 STD_SETTER_CMD_IMPL_F_S(XYInterpolationCurve, SetInterpolationData, XYInterpolationCurve::InterpolationData, interpolationData, recalculate);
94 void XYInterpolationCurve::setInterpolationData(const XYInterpolationCurve::InterpolationData& interpolationData) {
96  exec(new XYInterpolationCurveSetInterpolationDataCmd(d, interpolationData, ki18n("%1: set options and perform the interpolation")));
97 }
98 
99 //##############################################################################
100 //######################### Private implementation #############################
101 //##############################################################################
103 }
104 
105 //no need to delete xColumn and yColumn, they are deleted
106 //when the parent aspect is removed
108 
110  QElapsedTimer timer;
111  timer.start();
112 
113  //create interpolation result columns if not available yet, clear them otherwise
114  if (!xColumn) {
117  xVector = static_cast<QVector<double>* >(xColumn->data());
118  yVector = static_cast<QVector<double>* >(yColumn->data());
119 
120  xColumn->setHidden(true);
121  q->addChild(xColumn);
122  yColumn->setHidden(true);
123  q->addChild(yColumn);
124 
125  q->setUndoAware(false);
126  q->setXColumn(xColumn);
127  q->setYColumn(yColumn);
128  q->setUndoAware(true);
129  } else {
130  xVector->clear();
131  yVector->clear();
132  }
133 
134  // clear the previous result
136 
137  //determine the data source columns
138  const AbstractColumn* tmpXDataColumn = nullptr;
139  const AbstractColumn* tmpYDataColumn = nullptr;
141  //spreadsheet columns as data source
142  tmpXDataColumn = xDataColumn;
143  tmpYDataColumn = yDataColumn;
144  } else {
145  //curve columns as data source
146  tmpXDataColumn = dataSourceCurve->xColumn();
147  tmpYDataColumn = dataSourceCurve->yColumn();
148  }
149 
150  if (!tmpXDataColumn || !tmpYDataColumn) {
152  emit q->dataChanged();
154  return;
155  }
156 
157  //check column sizes
158  if (tmpXDataColumn->rowCount() != tmpYDataColumn->rowCount()) {
160  interpolationResult.valid = false;
161  interpolationResult.status = i18n("Number of x and y data points must be equal.");
163  emit q->dataChanged();
165  return;
166  }
167 
168  //copy all valid data point for the interpolation to temporary vectors
169  QVector<double> xdataVector;
170  QVector<double> ydataVector;
171 
172  double xmin;
173  double xmax;
175  xmin = tmpXDataColumn->minimum();
176  xmax = tmpXDataColumn->maximum();
177  } else {
178  xmin = interpolationData.xRange.first();
179  xmax = interpolationData.xRange.last();
180  }
181 
182  XYAnalysisCurve::copyData(xdataVector, ydataVector, tmpXDataColumn, tmpYDataColumn, xmin, xmax);
183 
184  //number of data points to interpolate
185  const size_t n = (size_t)xdataVector.size();
186  if (n < 2) {
188  interpolationResult.valid = false;
189  interpolationResult.status = i18n("Not enough data points available.");
191  emit q->dataChanged();
193  return;
194  }
195 
196  double* xdata = xdataVector.data();
197  double* ydata = ydataVector.data();
198 
199  // interpolation settings
202  const double tension = interpolationData.tension;
203  const double continuity = interpolationData.continuity;
204  const double bias = interpolationData.bias;
206  const size_t npoints = interpolationData.npoints;
207 
208  DEBUG("type:"<<nsl_interp_type_name[type]);
209  DEBUG("cubic Hermite variant:"<<nsl_interp_pch_variant_name[variant]<<tension<<continuity<<bias);
210  DEBUG("evaluate:"<<nsl_interp_evaluate_name[evaluate]);
211  DEBUG("npoints ="<<npoints);
212 
213 ///////////////////////////////////////////////////////////
214  int status = 0;
215 
216  gsl_interp_accel *acc = gsl_interp_accel_alloc();
217  gsl_spline *spline = nullptr;
218  switch (type) {
220  spline = gsl_spline_alloc(gsl_interp_linear, n);
221  status = gsl_spline_init(spline, xdata, ydata, n);
222  break;
224  spline = gsl_spline_alloc(gsl_interp_polynomial, n);
225  status = gsl_spline_init(spline, xdata, ydata, n);
226  break;
228  spline = gsl_spline_alloc(gsl_interp_cspline, n);
229  status = gsl_spline_init(spline, xdata, ydata, n);
230  break;
232  spline = gsl_spline_alloc(gsl_interp_cspline_periodic, n);
233  status = gsl_spline_init(spline, xdata, ydata, n);
234  break;
236  spline = gsl_spline_alloc(gsl_interp_akima, n);
237  status = gsl_spline_init(spline, xdata, ydata, n);
238  break;
240  spline = gsl_spline_alloc(gsl_interp_akima_periodic, n);
241  status = gsl_spline_init(spline, xdata, ydata, n);
242  break;
244 #if GSL_MAJOR_VERSION >= 2
245  spline = gsl_spline_alloc(gsl_interp_steffen, n);
246  status = gsl_spline_init(spline, xdata, ydata, n);
247 #endif
248  break;
250  case nsl_interp_type_pch:
253  break;
254  }
255 
256  xVector->resize((int)npoints);
257  yVector->resize((int)npoints);
258  for (unsigned int i = 0; i < npoints; i++) {
259  size_t a = 0, b = n-1;
260 
261  double x = xmin + i*(xmax-xmin)/(npoints-1);
262  (*xVector)[(int)i] = x;
263 
264  // find index a,b for interval [x[a],x[b]] around x[i] using bisection
266  while (b-a > 1) {
267  unsigned int j = floor((a+b)/2.);
268  if (xdata[j] > x)
269  b = j;
270  else
271  a = j;
272  }
273  }
274 
275  // evaluate interpolation
276  double t;
277  switch (type) {
285  switch (evaluate) {
287  (*yVector)[(int)i] = gsl_spline_eval(spline, x, acc);
288  break;
290  (*yVector)[(int)i] = gsl_spline_eval_deriv(spline, x, acc);
291  break;
293  (*yVector)[(int)i] = gsl_spline_eval_deriv2(spline, x, acc);
294  break;
296  (*yVector)[(int)i] = gsl_spline_eval_integ(spline, xmin, x, acc);
297  break;
298  }
299  break;
301  t = (x-xdata[a])/(xdata[b]-xdata[a]);
302  t = (1.-cos(M_PI*t))/2.;
303  (*yVector)[(int)i] = ydata[a] + t*(ydata[b]-ydata[a]);
304  break;
306  t = (x-xdata[a])/(xdata[b]-xdata[a]);
307  (*yVector)[(int)i] = ydata[a]*pow(ydata[b]/ydata[a],t);
308  break;
309  case nsl_interp_type_pch: {
310  t = (x-xdata[a])/(xdata[b]-xdata[a]);
311  double t2 = t*t, t3 = t2*t;
312  double h1 = 2.*t3-3.*t2+1, h2 = -2.*t3+3.*t2, h3 = t3-2*t2+t, h4 = t3-t2;
313  double m1 = 0.,m2 = 0.;
314  switch (variant) {
316  if (a == 0)
317  m1 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
318  else
319  m1 = ( (ydata[b]-ydata[a])/(xdata[b]-xdata[a]) + (ydata[a]-ydata[a-1])/(xdata[a]-xdata[a-1]) )/2.;
320  if (b == n-1)
321  m2 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
322  else
323  m2 = ( (ydata[b+1]-ydata[b])/(xdata[b+1]-xdata[b]) + (ydata[b]-ydata[a])/(xdata[b]-xdata[a]) )/2.;
324 
325  break;
327  if (a == 0)
328  m1 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
329  else
330  m1 = (ydata[b]-ydata[a-1])/(xdata[b]-xdata[a-1]);
331  if (b == n-1)
332  m2 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
333  else
334  m2 = (ydata[b+1]-ydata[a])/(xdata[b+1]-xdata[a]);
335 
336  break;
338  if (a == 0)
339  m1 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
340  else
341  m1 = (ydata[b]-ydata[a-1])/(xdata[b]-xdata[a-1]);
342  m1 *= (1.-tension);
343  if (b == n-1)
344  m2 = (ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
345  else
346  m2 = (ydata[b+1]-ydata[a])/(xdata[b+1]-xdata[a]);
347  m2 *= (1.-tension);
348 
349  break;
351  if (a == 0)
352  m1 = (1.+continuity)*(1.-bias)*(ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
353  else
354  m1 = ( (1.-continuity)*(1.+bias)*(ydata[a]-ydata[a-1])/(xdata[a]-xdata[a-1])
355  + (1.+continuity)*(1.-bias)*(ydata[b]-ydata[a])/(xdata[b]-xdata[a]) )/2.;
356  m1 *= (1.-tension);
357  if (b == n-1)
358  m2 = (1.+continuity)*(1.+bias)*(ydata[b]-ydata[a])/(xdata[b]-xdata[a]);
359  else
360  m2 = ( (1.+continuity)*(1.+bias)*(ydata[b]-ydata[a])/(xdata[b]-xdata[a])
361  + (1.-continuity)*(1.-bias)*(ydata[b+1]-ydata[b])/(xdata[b+1]-xdata[b]) )/2.;
362  m2 *= (1.-tension);
363 
364  break;
365  }
366 
367  // Hermite polynomial
368  (*yVector)[(int)i] = ydata[a]*h1+ydata[b]*h2+(xdata[b]-xdata[a])*(m1*h3+m2*h4);
369  }
370  break;
372  double v,dv;
373  nsl_interp_ratint(xdata, ydata, (int)n, x, &v, &dv);
374  (*yVector)[(int)i] = v;
375  //TODO: use error dv
376  break;
377  }
378  }
379  }
380 
381  // calculate "evaluate" option for own types
383  switch (evaluate) {
385  break;
387  nsl_diff_first_deriv_second_order(xVector->data(), yVector->data(), npoints);
388  break;
390  nsl_diff_second_deriv_second_order(xVector->data(), yVector->data(), npoints);
391  break;
393  nsl_int_trapezoid(xVector->data(), yVector->data(), npoints, 0);
394  break;
395  }
396  }
397 
398  // check values
399  for (int i = 0; i < (int)npoints; i++) {
400  if ((*yVector)[i] > std::numeric_limits<double>::max())
401  (*yVector)[i] = std::numeric_limits<double>::max();
402  else if ((*yVector)[i] < std::numeric_limits<double>::lowest())
403  (*yVector)[i] = std::numeric_limits<double>::lowest();
404  }
405 
406  gsl_spline_free(spline);
407  gsl_interp_accel_free(acc);
408 
409 ///////////////////////////////////////////////////////////
410 
411  //write the result
413  interpolationResult.valid = true;
414  interpolationResult.status = gslErrorToString(status);
415  interpolationResult.elapsedTime = timer.elapsed();
416 
417  //redraw the curve
419  emit q->dataChanged();
421 }
422 
423 //##############################################################################
424 //################## Serialization/Deserialization ###########################
425 //##############################################################################
426 //! Save as XML
427 void XYInterpolationCurve::save(QXmlStreamWriter* writer) const {
428  Q_D(const XYInterpolationCurve);
429 
430  writer->writeStartElement("xyInterpolationCurve");
431 
432  //write the base class
433  XYAnalysisCurve::save(writer);
434 
435  //write xy-interpolation-curve specific information
436  // interpolation data
437  writer->writeStartElement("interpolationData");
438  writer->writeAttribute( "autoRange", QString::number(d->interpolationData.autoRange) );
439  writer->writeAttribute( "xRangeMin", QString::number(d->interpolationData.xRange.first()) );
440  writer->writeAttribute( "xRangeMax", QString::number(d->interpolationData.xRange.last()) );
441  writer->writeAttribute( "type", QString::number(d->interpolationData.type) );
442  writer->writeAttribute( "variant", QString::number(d->interpolationData.variant) );
443  writer->writeAttribute( "tension", QString::number(d->interpolationData.tension) );
444  writer->writeAttribute( "continuity", QString::number(d->interpolationData.continuity) );
445  writer->writeAttribute( "bias", QString::number(d->interpolationData.bias) );
446  writer->writeAttribute( "npoints", QString::number(d->interpolationData.npoints) );
447  writer->writeAttribute( "pointsMode", QString::number(static_cast<int>(d->interpolationData.pointsMode)) );
448  writer->writeAttribute( "evaluate", QString::number(d->interpolationData.evaluate) );
449  writer->writeEndElement();// interpolationData
450 
451  // interpolation results (generated columns)
452  writer->writeStartElement("interpolationResult");
453  writer->writeAttribute( "available", QString::number(d->interpolationResult.available) );
454  writer->writeAttribute( "valid", QString::number(d->interpolationResult.valid) );
455  writer->writeAttribute( "status", d->interpolationResult.status );
456  writer->writeAttribute( "time", QString::number(d->interpolationResult.elapsedTime) );
457 
458  //save calculated columns if available
459  if (d->xColumn) {
460  d->xColumn->save(writer);
461  d->yColumn->save(writer);
462  }
463  writer->writeEndElement(); //"interpolationResult"
464 
465  writer->writeEndElement(); //"xyInterpolationCurve"
466 }
467 
468 //! Load from XML
469 bool XYInterpolationCurve::load(XmlStreamReader* reader, bool preview) {
471 
472  KLocalizedString attributeWarning = ki18n("Attribute '%1' missing or empty, default value is used");
473  QXmlStreamAttributes attribs;
474  QString str;
475 
476  while (!reader->atEnd()) {
477  reader->readNext();
478  if (reader->isEndElement() && reader->name() == "xyInterpolationCurve")
479  break;
480 
481  if (!reader->isStartElement())
482  continue;
483 
484  if (reader->name() == "xyAnalysisCurve") {
485  if ( !XYAnalysisCurve::load(reader, preview) )
486  return false;
487  } else if (!preview && reader->name() == "interpolationData") {
488  attribs = reader->attributes();
489  READ_INT_VALUE("autoRange", interpolationData.autoRange, bool);
490  READ_DOUBLE_VALUE("xRangeMin", interpolationData.xRange.first());
491  READ_DOUBLE_VALUE("xRangeMax", interpolationData.xRange.last());
492  READ_INT_VALUE("type", interpolationData.type, nsl_interp_type);
493  READ_INT_VALUE("variant", interpolationData.variant, nsl_interp_pch_variant);
494  READ_DOUBLE_VALUE("tension", interpolationData.tension);
495  READ_DOUBLE_VALUE("continuity", interpolationData.continuity);
496  READ_DOUBLE_VALUE("bias", interpolationData.bias);
497  READ_INT_VALUE("npoints", interpolationData.npoints, size_t);
498  READ_INT_VALUE("pointsMode", interpolationData.pointsMode, XYInterpolationCurve::PointsMode);
499  READ_INT_VALUE("evaluate", interpolationData.evaluate, nsl_interp_evaluate);
500  } else if (!preview && reader->name() == "interpolationResult") {
501  attribs = reader->attributes();
502  READ_INT_VALUE("available", interpolationResult.available, int);
506  } else if (reader->name() == "column") {
507  Column* column = new Column(QString(), AbstractColumn::ColumnMode::Numeric);
508  if (!column->load(reader, preview)) {
509  delete column;
510  return false;
511  }
512  if (column->name() == "x")
513  d->xColumn = column;
514  else if (column->name() == "y")
515  d->yColumn = column;
516  }
517  }
518 
519  if (preview)
520  return true;
521 
522  // wait for data to be read before using the pointers
523  QThreadPool::globalInstance()->waitForDone();
524 
525  if (d->xColumn && d->yColumn) {
526  d->xColumn->setHidden(true);
527  addChild(d->xColumn);
528 
529  d->yColumn->setHidden(true);
530  addChild(d->yColumn);
531 
532  d->xVector = static_cast<QVector<double>* >(d->xColumn->data());
533  d->yVector = static_cast<QVector<double>* >(d->yColumn->data());
534 
535  XYCurve::d_ptr->xColumn = d->xColumn;
536  XYCurve::d_ptr->yColumn = d->yColumn;
537 
539  }
540 
541  return true;
542 }
AspectType
STD_SETTER_CMD_IMPL_F_S(XYInterpolationCurve, SetInterpolationData, XYInterpolationCurve::InterpolationData, interpolationData, recalculate)
void setUndoAware(bool)
void addChild(AbstractAspect *)
Add the given Aspect to my list of children.
QString name() const
void exec(QUndoCommand *)
Execute the given command, pushing it on the undoStack() if available.
AbstractAspectPrivate * d
void setHidden(bool)
Set "hidden" property, i.e. whether to exclude this aspect from being shown in the explorer.
Interface definition for data with column logic.
virtual int rowCount() const =0
Return the data vector size.
virtual double maximum(int count=0) const
virtual double minimum(int count=0) const
Aspect that manages a column.
Definition: Column.h:42
void * data() const
Definition: Column.cpp:852
bool load(XmlStreamReader *, bool preview) override
Load the column from XML.
Definition: Column.cpp:1152
XYAnalysisCurve::DataSourceType dataSourceType
const AbstractColumn * yDataColumn
QVector< double > * xVector
const AbstractColumn * xDataColumn
QVector< double > * yVector
Base class for all analysis curves.
static void copyData(QVector< double > &xData, QVector< double > &yData, const AbstractColumn *xDataColumn, const AbstractColumn *yDataColumn, double xMin, double xMax)
void save(QXmlStreamWriter *) const override
Save as XML.
bool load(XmlStreamReader *, bool preview) override
Load from XML.
bool sourceDataChangedSinceLastRecalc
const AbstractColumn * xColumn
void recalcLogicalPoints()
Definition: XYCurve.cpp:1090
const AbstractColumn * yColumn
void dataChanged()
XYCurvePrivate *const d_ptr
Definition: XYCurve.h:186
void recalcLogicalPoints()
Definition: XYCurve.cpp:765
XYInterpolationCurve::InterpolationData interpolationData
XYInterpolationCurvePrivate(XYInterpolationCurve *)
XYInterpolationCurve::InterpolationResult interpolationResult
~XYInterpolationCurvePrivate() override
A xy-curve defined by an interpolation.
~XYInterpolationCurve() override
void save(QXmlStreamWriter *) const override
Save as XML.
bool load(XmlStreamReader *, bool preview) override
Load from XML.
QIcon icon() const override
const InterpolationResult & interpolationResult() const
XYInterpolationCurve(const QString &name)
XML stream parser that supports errors as well as warnings. This class also adds line and column numb...
#define READ_DOUBLE_VALUE(name, var)
Definition: macros.h:475
#define DEBUG(x)
Definition: macros.h:50
#define READ_INT_VALUE(name, var, type)
Definition: macros.h:468
#define BASIC_SHARED_D_READER_IMPL(classname, type, method, var)
Definition: macros.h:127
#define READ_STRING_VALUE(name, var)
Definition: macros.h:482
class Origin::Variant variant
#define i18n(m)
Definition: nsl_common.h:38
int nsl_diff_first_deriv_second_order(const double *x, double *y, const size_t n)
Definition: nsl_diff.c:82
int nsl_diff_second_deriv_second_order(const double *x, double *y, const size_t n)
Definition: nsl_diff.c:218
int nsl_int_trapezoid(const double *x, double *y, const size_t n, int abs)
Definition: nsl_int.c:59
int nsl_interp_ratint(double *x, double *y, int n, double xn, double *v, double *dv)
Definition: nsl_interp.c:37
const char * nsl_interp_pch_variant_name[]
Definition: nsl_interp.c:34
const char * nsl_interp_type_name[]
Definition: nsl_interp.c:31
const char * nsl_interp_evaluate_name[]
Definition: nsl_interp.c:35
nsl_interp_evaluate
Definition: nsl_interp.h:44
@ nsl_interp_evaluate_function
Definition: nsl_interp.h:44
@ nsl_interp_evaluate_derivative
Definition: nsl_interp.h:44
@ nsl_interp_evaluate_second_derivative
Definition: nsl_interp.h:44
@ nsl_interp_evaluate_integral
Definition: nsl_interp.h:45
nsl_interp_pch_variant
Definition: nsl_interp.h:39
@ nsl_interp_pch_variant_kochanek_bartels
Definition: nsl_interp.h:40
@ nsl_interp_pch_variant_finite_difference
Definition: nsl_interp.h:39
@ nsl_interp_pch_variant_catmull_rom
Definition: nsl_interp.h:39
@ nsl_interp_pch_variant_cardinal
Definition: nsl_interp.h:39
nsl_interp_type
Definition: nsl_interp.h:33
@ nsl_interp_type_akima_periodic
Definition: nsl_interp.h:34
@ nsl_interp_type_polynomial
Definition: nsl_interp.h:33
@ nsl_interp_type_akima
Definition: nsl_interp.h:34
@ nsl_interp_type_cspline_periodic
Definition: nsl_interp.h:33
@ nsl_interp_type_linear
Definition: nsl_interp.h:33
@ nsl_interp_type_rational
Definition: nsl_interp.h:35
@ nsl_interp_type_cosine
Definition: nsl_interp.h:34
@ nsl_interp_type_pch
Definition: nsl_interp.h:35
@ nsl_interp_type_cspline
Definition: nsl_interp.h:33
@ nsl_interp_type_exponential
Definition: nsl_interp.h:35
@ nsl_interp_type_steffen
Definition: nsl_interp.h:34
XYInterpolationCurve::PointsMode pointsMode