Skip to content

Commit

Permalink
Use only knots inside of support for spacing bound
Browse files Browse the repository at this point in the history
It is perfectly legitimate for the first and last `order` knots to
be identical.  Exclude these from the minimum and maximum knot
spacing calculation to prevent elements of rmin_sep from becoming
infinite.
  • Loading branch information
jvansanten committed Apr 15, 2024
1 parent 05f4b1c commit a0e5a24
Show file tree
Hide file tree
Showing 6 changed files with 11 additions and 11 deletions.
2 changes: 1 addition & 1 deletion include/photospline/detail/bspline_eval.h
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,7 @@ bool splinetable<Alloc>::searchcenters(const double* x, int* centers) const
}

uint32_t min = order[i];
uint32_t max = nknots[i]-2;
uint32_t max = naxes[i];
double diff = x[i] - knots[i][min];
uint32_t hi = diff*rmin_sep[i];
uint32_t lo = diff*rmax_sep[i];
Expand Down
2 changes: 1 addition & 1 deletion include/photospline/detail/convolve.h
Original file line number Diff line number Diff line change
Expand Up @@ -146,7 +146,7 @@ void splinetable<Alloc>::convolve(const uint32_t dim, const double* conv_knots,
}

// Update knot separations
dknot_bounds();
fill_knot_spacing_bounds();

/*
* NB: A monotonic function remains monotonic after convolution
Expand Down
2 changes: 1 addition & 1 deletion include/photospline/detail/fit.h
Original file line number Diff line number Diff line change
Expand Up @@ -104,7 +104,7 @@ void splinetable<Alloc>::fit(const ::ndsparse& data,
std::copy(knots[i].begin(),knots[i].end(),this->knots[i]);
dummy_knots[i]=&this->knots[i][0];
}
dknot_bounds();
fill_knot_spacing_bounds();
//same deal for the coordinates
std::unique_ptr<const double*[]> dummy_coords(new const double*[ndim]);
for(uint32_t i=0; i<ndim; i++)
Expand Down
2 changes: 1 addition & 1 deletion include/photospline/detail/fitsio.h
Original file line number Diff line number Diff line change
Expand Up @@ -382,7 +382,7 @@ bool splinetable<Alloc>::read_fits_core(fitsfile* fits, const std::string& fileP
}
}

dknot_bounds();
fill_knot_spacing_bounds();

if(error!=0)
throw std::runtime_error("Error reading "+filePath+": Error "+std::to_string(error));
Expand Down
2 changes: 1 addition & 1 deletion include/photospline/detail/permute.h
Original file line number Diff line number Diff line change
Expand Up @@ -84,7 +84,7 @@ void splinetable<Alloc>::permuteDimensions(const std::vector<size_t>& permutatio
std::copy(t_coefficients.get(),t_coefficients.get()+ncoeffs,coefficients);

// Update knot separations
dknot_bounds();
fill_knot_spacing_bounds();
}

} //namespace photospline
Expand Down
12 changes: 6 additions & 6 deletions include/photospline/splinetable.h
Original file line number Diff line number Diff line change
Expand Up @@ -277,16 +277,16 @@ class splinetable{
lastKnots[nknots[inputDim]-1]=2*lastKnots[nknots[inputDim]-2]-lastKnots[nknots[inputDim]-3];
}

rmin_sep=allocate<double>(ndim);
rmax_sep=allocate<double>(ndim);
dknot_bounds();

//set naxes
naxes=allocate<uint64_t>(ndim);
for(unsigned int i=0; i<inputDim; i++)
naxes[i] = tables.front()->get_ncoeffs(i);
naxes[inputDim]=tables.size();

rmin_sep=allocate<double>(ndim);
rmax_sep=allocate<double>(ndim);
fill_knot_spacing_bounds();

//copy coefficients
unsigned long nCoeffs=std::accumulate(naxes, naxes+ndim, 1UL, std::multiplies<uint64_t>());
unsigned long nInputCoeffs=std::accumulate(naxes, naxes+ndim-1, 1UL, std::multiplies<uint64_t>());
Expand Down Expand Up @@ -856,10 +856,10 @@ class splinetable{
///Write to a file
void write_fits_core(fitsfile*) const;

void dknot_bounds() {
void fill_knot_spacing_bounds() {
for (uint32_t i = 0; i < ndim; i++) {
uint32_t min = order[i];
uint32_t max = nknots[i]-2;
uint32_t max = naxes[i];
double mini = DBL_MAX;
double maxi = 0;
for (uint32_t j = min; j < max; j++) {
Expand Down

0 comments on commit a0e5a24

Please sign in to comment.