Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Prefactoring: Rename and move interfaces #341

Merged
merged 11 commits into from
Sep 18, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ ds = FakeDataSet(512, 512, 10)
imgs = ds.stack.get_images()

# Get the timestamp of the first image.
t0 = imgs[0].get_time()
t0 = imgs[0].get_obstime()
print(f"Image times start at {t0}.")

# Specify an artificial object
Expand All @@ -110,7 +110,7 @@ velocity = (2, 0)

# Inject object into images
for im in imgs:
dt = im.get_time() - t0
dt = im.get_obstime() - t0
im.add_object(position[0] + dt * velocity[0],
position[1] + dt * velocity[1],
flux)
Expand Down
8 changes: 4 additions & 4 deletions notebooks/Kbmod_Reference.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -167,7 +167,7 @@
"outputs": [],
"source": [
"im = kb.layered_image(im_path + \"000000.fits\", p)\n",
"print(f\"Loaded a {im.get_width()} by {im.get_height()} image at time {im.get_time()}\")"
"print(f\"Loaded a {im.get_width()} by {im.get_height()} image at time {im.get_obstime()}\")"
]
},
{
Expand Down Expand Up @@ -204,8 +204,8 @@
"source": [
"print(f\"Width = {im.get_width()}\")\n",
"print(f\"Height = {im.get_height()}\")\n",
"print(f\"Pixels Per Image = {im.get_ppi()}\")\n",
"print(f\"Time = {im.get_time()}\")"
"print(f\"Pixels Per Image = {im.get_npixels()}\")\n",
"print(f\"Time = {im.get_obstime()}\")"
]
},
{
Expand Down Expand Up @@ -461,7 +461,7 @@
"# where the stack contains images of different sizes.\n",
"print(f\"Width = {stack.get_width()}\")\n",
"print(f\"Height = {stack.get_height()}\")\n",
"print(f\"Pixels Per Image = {stack.get_ppi()}\")\n",
"print(f\"Pixels Per Image = {stack.get_npixels()}\")\n",
"\n",
"# Retrieve a list of layered_images back from the stack.\n",
"stack.get_images()\n",
Expand Down
4 changes: 2 additions & 2 deletions notebooks/create_fake_data.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -111,9 +111,9 @@
],
"source": [
"# Check the object was inserted correctly.\n",
"t0 = ds.stack.get_single_image(0).get_time()\n",
"t0 = ds.stack.get_single_image(0).get_obstime()\n",
"for i in range(ds.stack.img_count()):\n",
" ti = ds.stack.get_single_image(i).get_time()\n",
" ti = ds.stack.get_single_image(i).get_obstime()\n",
" dt = ti - t0\n",
" px = int(trj.x + dt * trj.x_v + 0.5)\n",
" py = int(trj.y + dt * trj.y_v + 0.5)\n",
Expand Down
2 changes: 1 addition & 1 deletion notebooks/kbmod_visualize.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,7 @@
"source": [
"p = kb.psf(1.0)\n",
"im = kb.layered_image(im_path + \"000000.fits\", p)\n",
"print(f\"Loaded a {im.get_width()} by {im.get_height()} image at time {im.get_time()}\")"
"print(f\"Loaded a {im.get_width()} by {im.get_height()} image at time {im.get_obstime()}\")"
]
},
{
Expand Down
2 changes: 1 addition & 1 deletion src/kbmod/analysis_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -130,7 +130,7 @@ def load_images(
if verbose:
print(f"Loading file: {full_file_path}")
img = kb.layered_image(full_file_path, psf)
img.set_time(time_stamp)
img.set_obstime(time_stamp)

# Save the file, time, and image information.
img_info.append(header_info)
Expand Down
2 changes: 1 addition & 1 deletion src/kbmod/image_info.py
Original file line number Diff line number Diff line change
Expand Up @@ -140,7 +140,7 @@ def set_epoch(self, epoch):
"""
self.epoch_ = Time(epoch)
if self.image is not None:
self.image.set_time(self.epoch_.mjd)
self.image.set_obstime(self.epoch_.mjd)
self.epoch_set_ = True

def get_epoch(self, none_if_unset=False):
Expand Down
12 changes: 6 additions & 6 deletions src/kbmod/search/ImageStack.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ void ImageStack::extractImageTimes() {
// Load image times
imageTimes = std::vector<float>();
for (auto& i : images) {
imageTimes.push_back(float(i.getTime()));
imageTimes.push_back(float(i.getObstime()));
}
}

Expand Down Expand Up @@ -85,7 +85,7 @@ void ImageStack::convolvePSF() {
for (auto& i : images) i.convolvePSF();
}

void ImageStack::saveGlobalMask(const std::string& path) { globalMask.saveToFile(path, false); }
void ImageStack::saveGlobalMask(const std::string& path) { globalMask.saveToFile(path); }

void ImageStack::saveImages(const std::string& path) {
for (auto& i : images) i.saveLayers(path);
Expand Down Expand Up @@ -115,21 +115,21 @@ void ImageStack::growMask(int steps) {
}

void ImageStack::createGlobalMask(int flags, int threshold) {
int ppi = getPPI();
int npixels = getNPixels();

// For each pixel count the number of images where it is masked.
std::vector<int> counts(ppi, 0);
std::vector<int> counts(npixels, 0);
for (unsigned int img = 0; img < images.size(); ++img) {
float* imgMask = images[img].getMDataRef();
// Count the number of times a pixel has any of the flags
for (unsigned int pixel = 0; pixel < ppi; ++pixel) {
for (unsigned int pixel = 0; pixel < npixels; ++pixel) {
if ((flags & static_cast<int>(imgMask[pixel])) != 0) counts[pixel]++;
}
}

// Set all pixels below threshold to 0 and all above to 1
float* globalM = globalMask.getDataRef();
for (unsigned int p = 0; p < ppi; ++p) {
for (unsigned int p = 0; p < npixels; ++p) {
globalM[p] = counts[p] < threshold ? 0.0 : 1.0;
}
}
Expand Down
2 changes: 1 addition & 1 deletion src/kbmod/search/ImageStack.h
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@ class ImageStack {
unsigned imgCount() const { return images.size(); }
unsigned getWidth() const { return images.size() > 0 ? images[0].getWidth() : 0; }
unsigned getHeight() const { return images.size() > 0 ? images[0].getHeight() : 0; }
unsigned getPPI() const { return images.size() > 0 ? images[0].getPPI() : 0; }
unsigned getNPixels() const { return images.size() > 0 ? images[0].getNPixels() : 0; }
std::vector<LayeredImage>& getImages() { return images; }
const std::vector<float>& getTimes() const { return imageTimes; }
float* getTimesDataRef() { return imageTimes.data(); }
Expand Down
10 changes: 5 additions & 5 deletions src/kbmod/search/KBMOSearch.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -223,8 +223,8 @@ void KBMOSearch::saveImages(const std::string& path) {
std::string number = std::to_string(i);
// Add leading zeros
number = std::string(4 - number.length(), '0') + number;
psiImages[i].saveToFile(path + "/psi/PSI" + number + ".fits", false);
phiImages[i].saveToFile(path + "/phi/PHI" + number + ".fits", false);
psiImages[i].saveToFile(path + "/psi/PSI" + number + ".fits");
phiImages[i].saveToFile(path + "/phi/PHI" + number + ".fits");
}
}

Expand Down Expand Up @@ -262,10 +262,10 @@ void KBMOSearch::fillPsiAndPhiVects(const std::vector<RawImage>& psiImgs,
assert(num_images > 0);
assert(phiImgs.size() == num_images);

int num_pixels = psiImgs[0].getPPI();
int num_pixels = psiImgs[0].getNPixels();
for (int i = 0; i < num_images; ++i) {
assert(psiImgs[i].getPPI() == num_pixels);
assert(phiImgs[i].getPPI() == num_pixels);
assert(psiImgs[i].getNPixels() == num_pixels);
assert(phiImgs[i].getNPixels() == num_pixels);
}

psiVect->clear();
Expand Down
115 changes: 33 additions & 82 deletions src/kbmod/search/LayeredImage.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -15,14 +15,25 @@ LayeredImage::LayeredImage(std::string path, const PointSpreadFunc& psf) : psf(p
int fBegin = path.find_last_of("/");
int fEnd = path.find_last_of(".fits") - 4;
fileName = path.substr(fBegin, fEnd - fBegin);
readHeader(path);
science = RawImage(width, height);
mask = RawImage(width, height);
variance = RawImage(width, height);
loadLayers(path);

science = RawImage();
science.loadFromFile(path, 1);
width = science.getWidth();
height = science.getHeight();

mask = RawImage();
mask.loadFromFile(path, 2);

variance = RawImage();
variance.loadFromFile(path, 3);

if (width != variance.getWidth() or height != variance.getHeight())
throw std::runtime_error("Science and Variance layers are not the same size.");
if (width != mask.getWidth() or height != mask.getHeight())
throw std::runtime_error("Science and Mask layers are not the same size.");
}

LayeredImage::LayeredImage(const RawImage& sci, const RawImage& var, const RawImage& msk, float time,
LayeredImage::LayeredImage(const RawImage& sci, const RawImage& var, const RawImage& msk,
const PointSpreadFunc& psf)
: psf(psf), psfSQ(psf) {
// Get the dimensions of the science layer and check for consistency with
Expand All @@ -35,7 +46,6 @@ LayeredImage::LayeredImage(const RawImage& sci, const RawImage& var, const RawIm
throw std::runtime_error("Science and Mask layers are not the same size.");

// Set the remaining variables.
captureTime = time;
psfSQ.squarePSF();

// Copy the image layers.
Expand All @@ -54,7 +64,6 @@ LayeredImage::LayeredImage(std::string name, int w, int h, float noiseStDev, flo
fileName = name;
width = w;
height = h;
captureTime = time;
psfSQ.squarePSF();

std::vector<float> rawSci(width * height);
Expand All @@ -65,62 +74,14 @@ LayeredImage::LayeredImage(std::string name, int w, int h, float noiseStDev, flo
}
std::normal_distribution<float> distrib(0.0, noiseStDev);
for (float& p : rawSci) p = distrib(generator);

science = RawImage(w, h, rawSci);
science.setObstime(time);

mask = RawImage(w, h, std::vector<float>(w * h, 0.0));
variance = RawImage(w, h, std::vector<float>(w * h, pixelVariance));
}

/* Read the image dimensions and capture time from header */
void LayeredImage::readHeader(const std::string& filePath) {
fitsfile* fptr;
int status = 0;
int mjdStatus = 0;
int fileNotFound;

// Open header to read MJD
if (fits_open_file(&fptr, filePath.c_str(), READONLY, &status))
throw std::runtime_error("Could not open file");

// Read image capture time, ignore error if does not exist
captureTime = 0.0;
fits_read_key(fptr, TDOUBLE, "MJD", &captureTime, NULL, &mjdStatus);

if (fits_close_file(fptr, &status)) fits_report_error(stderr, status);

// Reopen header for first layer to get image dimensions
if (fits_open_file(&fptr, (filePath + "[1]").c_str(), READONLY, &status))
fits_report_error(stderr, status);

// Read image Dimensions
long dimensions[2];
if (fits_read_keys_lng(fptr, "NAXIS", 1, 2, dimensions, &fileNotFound, &status))
fits_report_error(stderr, status);

width = dimensions[0];
height = dimensions[1];

if (fits_close_file(fptr, &status)) fits_report_error(stderr, status);
}

void LayeredImage::loadLayers(const std::string& filePath) {
// Load images from file into layers' pixels
readFitsImg((filePath + "[1]").c_str(), science.getDataRef());
readFitsImg((filePath + "[2]").c_str(), mask.getDataRef());
readFitsImg((filePath + "[3]").c_str(), variance.getDataRef());
}

void LayeredImage::readFitsImg(const char* name, float* target) {
fitsfile* fptr;
int nullval = 0;
int anynull;
int status = 0;

if (fits_open_file(&fptr, name, READONLY, &status)) fits_report_error(stderr, status);
if (fits_read_img(fptr, TFLOAT, 1, getPPI(), &nullval, target, &anynull, &status))
fits_report_error(stderr, status);
if (fits_close_file(fptr, &status)) fits_report_error(stderr, status);
}

void LayeredImage::setPSF(const PointSpreadFunc& new_psf) {
psf = new_psf;
psfSQ = new_psf;
Expand Down Expand Up @@ -165,7 +126,7 @@ void LayeredImage::applyGlobalMask(const RawImage& globalM) {
}

void LayeredImage::applyMaskThreshold(float thresh) {
const int numPixels = getPPI();
const int numPixels = getNPixels();
float* sciPix = science.getDataRef();
float* varPix = variance.getDataRef();
for (int i = 0; i < numPixels; ++i) {
Expand All @@ -178,13 +139,13 @@ void LayeredImage::applyMaskThreshold(float thresh) {

void LayeredImage::subtractTemplate(const RawImage& subTemplate) {
assert(getHeight() == subTemplate.getHeight() && getWidth() == subTemplate.getWidth());
const int numPixels = getPPI();
const int numPixels = getNPixels();

float* sciPix = science.getDataRef();
const std::vector<float>& tempPix = subTemplate.getPixels();
const std::vector<float>& temNPixelsx = subTemplate.getPixels();
for (unsigned i = 0; i < numPixels; ++i) {
if ((sciPix[i] != NO_DATA) && (tempPix[i] != NO_DATA)) {
sciPix[i] -= tempPix[i];
if ((sciPix[i] != NO_DATA) && (temNPixelsx[i] != NO_DATA)) {
sciPix[i] -= temNPixelsx[i];
}
}
}
Expand All @@ -193,6 +154,8 @@ void LayeredImage::saveLayers(const std::string& path) {
fitsfile* fptr;
int status = 0;
long naxes[2] = {0, 0};
double obstime = science.getObstime();

fits_create_file(&fptr, (path + fileName + ".fits").c_str(), &status);

// If we are unable to create the file, check if it already exists
Expand All @@ -207,25 +170,13 @@ void LayeredImage::saveLayers(const std::string& path) {
}

fits_create_img(fptr, SHORT_IMG, 0, naxes, &status);
fits_update_key(fptr, TDOUBLE, "MJD", &captureTime, "[d] Generated Image time", &status);
fits_update_key(fptr, TDOUBLE, "MJD", &obstime, "[d] Generated Image time", &status);
fits_close_file(fptr, &status);
fits_report_error(stderr, status);

science.saveToFile(path + fileName + ".fits", true);
mask.saveToFile(path + fileName + ".fits", true);
variance.saveToFile(path + fileName + ".fits", true);
}

void LayeredImage::saveSci(const std::string& path) {
science.saveToFile(path + fileName + "SCI.fits", false);
}

void LayeredImage::saveMask(const std::string& path) {
mask.saveToFile(path + fileName + "MASK.fits", false);
}

void LayeredImage::saveVar(const std::string& path) {
variance.saveToFile(path + fileName + "VAR.fits", false);
science.appendLayerToFile(path + fileName + ".fits");
mask.appendLayerToFile(path + fileName + ".fits");
variance.appendLayerToFile(path + fileName + ".fits");
}

void LayeredImage::setScience(RawImage& im) {
Expand Down Expand Up @@ -255,7 +206,7 @@ RawImage LayeredImage::generatePsiImage() {
float* varArray = getVDataRef();

// Set each of the result pixels.
const int num_pixels = getPPI();
const int num_pixels = getNPixels();
for (int p = 0; p < num_pixels; ++p) {
float varPix = varArray[p];
if (varPix != NO_DATA) {
Expand All @@ -277,7 +228,7 @@ RawImage LayeredImage::generatePhiImage() {
float* varArray = getVDataRef();

// Set each of the result pixels.
const int num_pixels = getPPI();
const int num_pixels = getNPixels();
for (int p = 0; p < num_pixels; ++p) {
float varPix = varArray[p];
if (varPix != NO_DATA) {
Expand Down
Loading