68 if(start<sp.
start){
return 1;}
75 else if (d>0) {
return 1;}
97 if(polyCount){free(polys);}
106 for( i=0 ; i<count ; i++ )
108 if( !c || sps[i].start!=polys[c-1].start ) polys[c++] = sps[i];
109 else{polys[c-1].
p+=sps[i].
p;}
113 template <
int Degree>
119 if(polyCount){free(polys);}
143 template<
int Degree2>
146 for(
int i=0;i<int(polyCount);i++){
147 polys[i].start=p.
polys[i].start;
148 polys[i].p=p.
polys[i].p;
157 for(
int i=0 ; i<int(polyCount) && t>polys[i].start ; i++ ) v+=polys[i].p(t);
165 double start,end,s,v=0;
173 for(
int i=0;i<int(polyCount) && polys[i].start<end;i++){
174 if(start<polys[i].start){s=polys[i].start;}
176 v+=polys[i].p.integral(s,end);
192 else if (i>=
int( polyCount)-1) {q.
polys[idx]=p.
polys[++j];}
193 else if(polys[i+1].start<p.
polys[j+1].start){q.
polys[idx]= polys[++i];}
209 else if (i>=
int( polyCount)-1) {q.
polys[idx].start=p.
polys[++j].start;q.
polys[idx].p=p.
polys[j].p*(-1.0);}
210 else if(polys[i+1].start<p.
polys[j+1].start){q.
polys[idx]= polys[++i];}
220 size_t idx=0,cnt=0,oldPolyCount=polyCount;
225 while(cnt<polyCount){
226 if (j>=
int( p.
polyCount)-1) {polys[idx]=oldPolys[++i];}
227 else if (i>=
int(oldPolyCount)-1) {polys[idx].
start= p.
polys[++j].start;polys[idx].p=p.
polys[j].p*scale;}
228 else if (oldPolys[i+1].start<p.
polys[j+1].start){polys[idx]=oldPolys[++i];}
229 else {polys[idx].
start= p.
polys[++j].start;polys[idx].p=p.
polys[j].p*scale;}
230 if(idx && polys[idx].start==polys[idx-1].start) {polys[idx-1].p+=polys[idx].p;}
239 template<
int Degree2>
243 int i,j,spCount=int(polyCount*p.
polyCount);
246 for(i=0;i<int(polyCount);i++){
256 template<
int Degree2>
260 for(
int i=0;i<int(polyCount);i++){
261 q.
polys[i].start=polys[i].start;
262 q.
polys[i].p=polys[i].p*p;
271 for(
size_t i=0;i<polyCount;i++){q.
polys[i]=polys[i].scale(s);}
279 for(
size_t i=0;i<polyCount;i++){q.
polys[i]=polys[i].shift(s);}
286 for(
size_t i=0;i<polyCount;i++){
287 q.polys[i].start=polys[i].start;
288 q.polys[i].p=polys[i].p.derivative();
297 for(i=0;i<int(polyCount);i++){
298 q.
polys[i].start=polys[i].start;
299 q.
polys[i].p=polys[i].p.integral();
311 for(
int i=0;i<int(polyCount);i++){polys[i].p*=s;}
317 for(
size_t i=0;i<polyCount;i++){polys[i].p/=s;}
355 printf(
"[-Infinity,Infinity]\n");
358 for(
size_t i=0;i<polyCount;i++){
360 if (polys[i ].start== DBL_MAX){printf(
"Infinity,");}
361 else if (polys[i ].start==-DBL_MAX){printf(
"-Infinity,");}
362 else {printf(
"%f,",polys[i].start);}
363 if(i+1==polyCount) {printf(
"Infinity]\t");}
364 else if (polys[i+1].start== DBL_MAX){printf(
"Infinity]\t");}
365 else if (polys[i+1].start==-DBL_MAX){printf(
"-Infinity]\t");}
366 else {printf(
"%f]\t",polys[i+1].start);}
379 q.
polys[0].start=-radius;
380 q.
polys[1].start= radius;
382 q.
polys[0].p.coefficients[0]= 1.0;
383 q.
polys[1].p.coefficients[0]=-1.0;
386 template<
int Degree >
400 for(
int i=0;i<int(polyCount);i++){
401 sps[2*i ].
start=polys[i].start-radius;
402 sps[2*i+1].
start=polys[i].start+radius;
403 p=polys[i].p.
integral()-polys[i].p.integral()(polys[i].start);
404 sps[2*i ].
p=p.
shift(-radius);
405 sps[2*i+1].
p=p.
shift( radius)*-1;
407 A.
set(sps,
int(polyCount*2));
409 return A*1.0/(2*radius);
414 std::vector<double> tempRoots;
417 for(
size_t i=0;i<polyCount;i++){
419 if(polys[i].start>max){
break;}
420 if(i<polyCount-1 && polys[i+1].start<min){
continue;}
422 for(
size_t j=0;j<tempRoots.size();j++){
423 if(tempRoots[j]>polys[i].start && (i+1==polyCount || tempRoots[j]<=polys[i+1].start)){
424 if(tempRoots[j]>min && tempRoots[j]<max){roots.push_back(tempRoots[j]);}
432 fwrite(&samples,
sizeof(
int),1,fp);
433 for(
int i=0;i<samples;i++){
434 double x=min+i*(max-min)/(samples-1);
436 fwrite(&v,
sizeof(
float),1,fp);
PPolynomial operator+(const PPolynomial &p) const
PPolynomial & operator*=(double s)
double Integral(void) const
static PPolynomial BSpline(double radius=0.5)
PPolynomial operator-(const PPolynomial &p) const
double operator()(double t) const
This file defines compatibility wrappers for low level I/O functions.
int operator<(const StartingPolynomial &sp) const
void write(FILE *fp, int samples, double min, double max) const
PPolynomial & operator-=(double s)
PPolynomial< Degree+1 > MovingAverage(double radius)
PPolynomial & operator=(const PPolynomial &p)
StartingPolynomial< Degree+Degree2 > operator*(const StartingPolynomial< Degree2 > &p) const
PPolynomial & operator+=(double s)
PPolynomial scale(double s) const
PPolynomial & operator/=(double s)
PPolynomial< Degree+1 > integral(void) const
StartingPolynomial scale(double s) const
double integral(double tMin, double tMax) const
PPolynomial operator/(double s) const
PPolynomial< Degree+Degree2 > operator*(const Polynomial< Degree2 > &p) const
StartingPolynomial shift(double t) const
PPolynomial< Degree-1 > derivative(void) const
void getSolutions(double c, std::vector< double > &roots, double EPS) const
PPolynomial & addScaled(const PPolynomial &poly, double scale)
StartingPolynomial< Degree > * polys
void reset(size_t newSize)
void getSolutions(double c, std::vector< double > &roots, double EPS, double min=-DBL_MAX, double max=DBL_MAX) const
PPolynomial shift(double t) const
Polynomial shift(double t) const
static int Compare(const void *v1, const void *v2)