Stuermer / EchelleSimulator

GNU General Public License v3.0
11 stars 2 forks source link

--custom spectra not currently functional? #21

Closed abgibbs closed 3 years ago

abgibbs commented 3 years ago

I'm trying to input a custom, continuous spectrum using the --custom1 and --custom2 options, but can't get it to work. Is this a currently functional feature or was it a planned feature? The only documentation I see on it is in the --help.

I've tried this with both docker and a local installation and the results are the same. No matter what I input, using --custom1 or custom2 always leads me to either a segmentation fault or CCFits can't read. I don't know what the expected format is for the custom spectra, but based on the code, I'm assuming it is meant to have the same format as the phoenix files.

When I look at how the custom1 and 2 arguments call the CustomSpectra class, to me the variables don't seem to line up with what is listed with --help, which leads me to believe this is not yet a functional feature? If not, is there any suggestion on how to hack it?

Thanks in advance!

Stuermer commented 3 years ago

In fact it is implemented - but not documented at all. I apologize! After some trial and error: the correct format for the fits files is: a spectrum file (in - I believe - erg/s/cm^2/A ) and a separate wavelength file (in angstroms). Here is a minimal python script to generate a 'step' spectrum from 4000 - 10000 angstroms

import astropy.io.fits as fits
import numpy as np

# spectrum in erg/s/cm^2/A 
hdu = fits.PrimaryHDU(data=np.array(np.arange(100000)>25000, dtype=float))
hdu.writeto('test_spectrum.fits',overwrite=True)

# wavelength in Angstroms
hdu = fits.PrimaryHDU(data=np.linspace(4000, 10000, 100000))
hdu.writeto('test_wl.fits',overwrite=True)

To use it with Echelle++ you have to call it like this:

echellesimulator -s MaroonX -t 0.1 --custom2 path/to/test_spectrum.fits,/path/to/test_wl.fits,10,10

It's important, that there is no whitespace between any of the arguments of --custom2 The numbers at the end are V magnitude just like for the PHOENIX spectra. Don't ask me why you need to type it twice. Looking at the code it looks like a bug to me :-)

one comment: I stopped working on this package because it is painful to maintain and distribute due to the dependencies to cfitsio and hdf. I rewrote the package in python here. I encourage you to use this. It's not entirely complete yet, but my hope is also by switching to python that this makes it easier for others to contribute. It might be a bit slower. But a lot easier to implement new features.

one question: out of curiosity - on what OS are you running this and have you been able to compile it locally? I'm still struggling on my new machine/new OS to get the 'Release' version running - while the 'Debug' version seems to work fine.

abgibbs commented 3 years ago

Thanks for the quick response! I wasn't aware of the python version, but I probably will switch to that for ease. Otherwise, I'll try out your example. I currently have it compiled on Ubuntu 18.04. I had some issues compiling as well on my OS, but after messing with the dependencies for a while I recompiled and it worked. Unfortunately, I'm not sure what I did that actually made the difference.

Stuermer commented 3 years ago

The python version isn't 'released' yet, but it is in principle working already. That's why I also haven't written anything here in the README yet. Just be aware that a) I do work on this at the moment. I will try to have always a working version in the master branch. I can't promise that it will always be 'backwards' compatible, since I'm still in the process of refactoring code b) there is a slight addition to the .hdf models required compared to Echelle++, since I'm now supporting different shapes/distributions at the input slit. Therefore, there is another attribute 'field_shape' required for each 'fiber' within the .hdf model file. You can convert your .hdf model using this snippet

with h5py.File('/path/to/hdfmodel.hdf', 'a') as hf:
    for k in hf.keys():
        if 'fiber' in k:
            hf[k].attrs['field_shape'] = 'circular'
    hf.flush()

You can replace 'circular' with 'rectangular' or 'octagonal' or 'hexagonal' depending on what you want to do.

abgibbs commented 3 years ago

It seems like there may be a bug with custom spectra in the master branch, or possibly it is just on my end. I couldn't get the snippets of code you sent to run successfully. After generating exactly the same test_spectrum.fits and test_wl.fits, then using your call, I run into a segmentation fault. The same thing happens for docker, and all other commands I've tried seem to be working, which leads me to believe it is not something on my end.

Using gdb debugger I get: `Spectrograph file found: ../data/spectrographs/MaroonX.hdf Read in 33 Orders

Iterating over elements in the file Simulating custom spectra with magnitude = 10

Program received signal SIGSEGV, Segmentation fault. 0x0000555555637c42 in MatrixSimulator::set_wavelength(std::vector<double, std::allocator >) () `

Is the 'Debug' version you referenced available? I may try to debug the code myself as well.

Stuermer commented 3 years ago

hm - I'm trying to reproduce that error. Can you tell me what operating system you're using, what cmake version and what gcc ? I have a bunch of testing environments now that I can try.

There are 4 different build types, 'Debug' , 'Release', 'RelNoParallel' and 'RelWithDebInfo' (see CMakeList.txt). While 'Release' will deliver the fastest program for debugging I recommend to build echelle++ like this:

cmake ../ -DCMAKE_BUILD_TYPE=Debug

abgibbs commented 3 years ago

The versions: Ubuntu 18.04.5 LTS gcc (Ubuntu 7.5.0-3ubuntu1~18.04) 7.5.0 cmake version 3.10.2

Thanks for looking into this.

Stuermer commented 3 years ago

You're lucky. 18.04 is the the version where I get it running without issues. I updated installation instructions from a fresh 18.04 system. Also, there was a dependency missing in CMakeList. And finally, there was indeed a tiny bug in source.cpp which caused your above segmentation fault. Please pull the latest version and try again. I'll leave the issue open in the meantime...

abgibbs commented 3 years ago

Yep, that fixed it! From a brief look at a few different custom spectra, the output seems to be as expected. I'll close the issue now.