Line data Source code
1 : /* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 : /*
3 : * This file is part of the LibreOffice project.
4 : *
5 : * This Source Code Form is subject to the terms of the Mozilla Public
6 : * License, v. 2.0. If a copy of the MPL was not distributed with this
7 : * file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 : *
9 : * This file incorporates work covered by the following license notice:
10 : *
11 : * Licensed to the Apache Software Foundation (ASF) under one or more
12 : * contributor license agreements. See the NOTICE file distributed
13 : * with this work for additional information regarding copyright
14 : * ownership. The ASF licenses this file to you under the Apache
15 : * License, Version 2.0 (the "License"); you may not use this file
16 : * except in compliance with the License. You may obtain a copy of
17 : * the License at http://www.apache.org/licenses/LICENSE-2.0 .
18 : */
19 :
20 : #include "PolynomialRegressionCurveCalculator.hxx"
21 : #include "macros.hxx"
22 : #include "RegressionCalculationHelper.hxx"
23 :
24 : #include <cmath>
25 : #include <rtl/math.hxx>
26 : #include <rtl/ustrbuf.hxx>
27 : #include "gauss.hxx"
28 :
29 : using namespace com::sun::star;
30 :
31 :
32 : namespace chart
33 : {
34 :
35 0 : PolynomialRegressionCurveCalculator::PolynomialRegressionCurveCalculator()
36 0 : {}
37 :
38 0 : PolynomialRegressionCurveCalculator::~PolynomialRegressionCurveCalculator()
39 0 : {}
40 :
41 : // ____ XRegressionCurveCalculator ____
42 0 : void SAL_CALL PolynomialRegressionCurveCalculator::recalculateRegression(
43 : const uno::Sequence< double >& aXValues,
44 : const uno::Sequence< double >& aYValues )
45 : throw (uno::RuntimeException)
46 : {
47 0 : rtl::math::setNan(&m_fCorrelationCoeffitient);
48 :
49 : RegressionCalculationHelper::tDoubleVectorPair aValues(
50 0 : RegressionCalculationHelper::cleanup( aXValues, aYValues, RegressionCalculationHelper::isValid()));
51 :
52 0 : sal_Int32 aNoElements = mForceIntercept ? mDegree : mDegree + 1;
53 0 : sal_Int32 aNumberOfPowers = 2 * aNoElements - 1;
54 :
55 0 : std::vector<double> aPowers;
56 0 : aPowers.resize(aNumberOfPowers, 0.0);
57 :
58 0 : sal_Int32 aNoColumns = aNoElements;
59 0 : sal_Int32 aNoRows = aNoElements + 1;
60 :
61 0 : std::vector<double> aMatrix;
62 0 : aMatrix.resize(aNoColumns * aNoRows, 0.0);
63 :
64 0 : const size_t aNoValues = aValues.first.size();
65 :
66 0 : double yAverage = 0.0;
67 :
68 0 : for( size_t i = 0; i < aNoValues; ++i )
69 : {
70 0 : double x = aValues.first[i];
71 0 : double y = aValues.second[i];
72 :
73 0 : for (sal_Int32 j = 0; j < aNumberOfPowers; j++)
74 : {
75 0 : if (mForceIntercept)
76 0 : aPowers[j] += std::pow(x, (int) j + 2);
77 : else
78 0 : aPowers[j] += std::pow(x, (int) j);
79 : }
80 :
81 0 : for (sal_Int32 j = 0; j < aNoElements; j++)
82 : {
83 0 : if (mForceIntercept)
84 0 : aMatrix[j * aNoRows + aNoElements] += std::pow(x, (int) j + 1) * ( y - mInterceptValue );
85 : else
86 0 : aMatrix[j * aNoRows + aNoElements] += std::pow(x, (int) j) * y;
87 : }
88 :
89 0 : yAverage += y;
90 : }
91 :
92 0 : yAverage = yAverage / aNoValues;
93 :
94 0 : for (sal_Int32 y = 0; y < aNoElements; y++)
95 : {
96 0 : for (sal_Int32 x = 0; x < aNoElements; x++)
97 : {
98 0 : aMatrix[y * aNoRows + x] = aPowers[y + x];
99 : }
100 : }
101 :
102 0 : mResult.clear();
103 0 : mResult.resize(aNoElements, 0.0);
104 :
105 0 : solve(aMatrix, aNoColumns, aNoRows, mResult, 1.0e-20);
106 :
107 : // Set intercept value if force intercept is enabled
108 0 : if (mForceIntercept) {
109 0 : mResult.insert( mResult.begin(), mInterceptValue );
110 : }
111 :
112 : // Calculate correlation coeffitient
113 0 : double aSumError = 0.0;
114 0 : double aSumTotal = 0.0;
115 :
116 0 : for( size_t i = 0; i < aNoValues; ++i )
117 : {
118 0 : double x = aValues.first[i];
119 0 : double yActual = aValues.second[i];
120 0 : double yPredicted = getCurveValue( x );
121 0 : aSumTotal += (yActual - yAverage) * (yActual - yAverage);
122 0 : aSumError += (yActual - yPredicted) * (yActual - yPredicted);
123 : }
124 :
125 0 : double aRSquared = 1.0 - (aSumError / aSumTotal);
126 0 : if (aRSquared > 0.0)
127 0 : m_fCorrelationCoeffitient = std::sqrt(aRSquared);
128 : else
129 0 : m_fCorrelationCoeffitient = 0.0;
130 0 : }
131 :
132 0 : double SAL_CALL PolynomialRegressionCurveCalculator::getCurveValue( double x )
133 : throw (lang::IllegalArgumentException,
134 : uno::RuntimeException)
135 : {
136 : double fResult;
137 0 : rtl::math::setNan(&fResult);
138 :
139 0 : if (mResult.empty())
140 : {
141 0 : return fResult;
142 : }
143 :
144 0 : fResult = 0.0;
145 0 : for (size_t i = 0; i<mResult.size(); i++)
146 : {
147 0 : fResult += mResult[i] * std::pow(x, (int) i);
148 : }
149 0 : return fResult;
150 : }
151 :
152 0 : uno::Sequence< geometry::RealPoint2D > SAL_CALL PolynomialRegressionCurveCalculator::getCurveValues(
153 : double min, double max, sal_Int32 nPointCount,
154 : const uno::Reference< chart2::XScaling >& xScalingX,
155 : const uno::Reference< chart2::XScaling >& xScalingY,
156 : sal_Bool bMaySkipPointsInCalculation )
157 : throw (lang::IllegalArgumentException,
158 : uno::RuntimeException)
159 : {
160 :
161 0 : return RegressionCurveCalculator::getCurveValues( min, max, nPointCount, xScalingX, xScalingY, bMaySkipPointsInCalculation );
162 : }
163 :
164 0 : OUString PolynomialRegressionCurveCalculator::ImplGetRepresentation(
165 : const uno::Reference< util::XNumberFormatter >& xNumFormatter,
166 : sal_Int32 nNumberFormatKey ) const
167 : {
168 0 : OUStringBuffer aBuf( "f(x) = ");
169 :
170 0 : sal_Int32 aLastIndex = mResult.size() - 1;
171 0 : for (sal_Int32 i = aLastIndex; i >= 0; i--)
172 : {
173 0 : double aValue = mResult[i];
174 0 : if (aValue == 0.0)
175 : {
176 0 : continue;
177 : }
178 0 : else if (aValue < 0.0)
179 : {
180 0 : aBuf.appendAscii( " - " );
181 : }
182 : else
183 : {
184 0 : if (i != aLastIndex)
185 0 : aBuf.appendAscii( " + " );
186 : }
187 :
188 0 : aBuf.append( getFormattedString( xNumFormatter, nNumberFormatKey, std::abs( aValue ) ) );
189 :
190 0 : if(i > 0)
191 : {
192 0 : if (i == 1)
193 : {
194 0 : aBuf.appendAscii( "x" );
195 : }
196 : else
197 : {
198 0 : aBuf.appendAscii( "x^" );
199 0 : aBuf.append(i);
200 : }
201 : }
202 : }
203 :
204 0 : return aBuf.makeStringAndClear();
205 : }
206 :
207 : } // namespace chart
208 :
209 : /* vim:set shiftwidth=4 softtabstop=4 expandtab: */
|
