Issue specifying custom atmophere

[Nicholaswogan]

Hi. I'm just beginning to play around with Picaso. I wanted to compute a reflected light spectrum of Earth. So, following the "Getting Started : Basic Inputs and Outputs" in the docs, I put together a file for P, T and mixing ratios for Earth. I've attached it here: earth.pt.zip

I then run

opacity = jdi.opannection(wave_range=[0.3,2]) #lets just use all defaults
start_case = jdi.inputs()

#phase angle
start_case.phase_angle(0) #radians

#define gravity
start_case.gravity(gravity=9.81, gravity_unit=u.Unit('m/(s**2)')) #any astropy units available

#define star
start_case.star(opacity, 5000,0,4.0) #opacity db, pysynphot database, temp, metallicity,

start_case.atmosphere(filename='earth.pt', delim_whitespace=True)

df = start_case.spectrum(opacity)

I get the following error

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
Input In [8], in <cell line: 1>()
----> 1 df = start_case.spectrum(opacity)

File ~/miniconda3/envs/picaso/lib/python3.9/site-packages/picaso/justdoit.py:1980, in inputs.spectrum(self, opacityclass, calculation, dimension, full_output, plot_opacity, as_dict)
   1976     self.inputs['surface_reflect'] = 0 
   1977     self.inputs['hard_surface'] = 0 
-> 1980 return picaso(self, opacityclass,dimension=dimension,calculation=calculation,
   1981     full_output=full_output, plot_opacity=plot_opacity, as_dict=as_dict)

File ~/miniconda3/envs/picaso/lib/python3.9/site-packages/picaso/justdoit.py:139, in picaso(bundle, opacityclass, dimension, calculation, full_output, plot_opacity, as_dict)
    137 atm.get_mmw()
    138 atm.get_density()
--> 139 atm.get_altitude(p_reference = p_reference)#will calculate altitude if r and m are given (opposed to just g)
    140 atm.get_column_density()
    142 #gets both continuum and needed rayleigh cross sections 
    143 #relies on continuum molecules are added into the opacity 
    144 #database. Rayleigh molecules are all in `rayleigh.py` 

File ~/miniconda3/envs/picaso/lib/python3.9/site-packages/picaso/atmsetup.py:371, in ATMSETUP.get_altitude(self, p_reference, constant_gravity)
    368         gravity[i] = self.c.G * self.planet.mass / ( z[i] )**2  
    370     scale_h = self.c.k_b * tlevel[i] / (mmw[i] * gravity[i])
--> 371     dz[i] = scale_h*(np.log(plevel[i+1]/ plevel[i]))
    372     z[i-1] = z[i] + dz[i]
    374 self.level['z'] = z

IndexError: index 200 is out of bounds for axis 0 with size 200

Wondering if you could help spot the issue with my input file.

[natashabatalha]

You can easily fix this by changing the reference pressure. See the tutorial for Terrestrial Planets.

It's just a quick addition of a line:

opacity = jdi.opannection(wave_range=[0.3,2]) #lets just use all defaults
start_case = jdi.inputs()

#phase angle
start_case.phase_angle(0) #radians

#define gravity
start_case.gravity(gravity=9.81, gravity_unit=u.Unit('m/(s**2)')) #any astropy units available

#define star
start_case.star(opacity, 5000,0,4.0) #opacity db, pysynphot database, temp, metallicity,

start_case.atmosphere(filename='earth.pt', delim_whitespace=True)

#NEW ADDITION
start_case.inputs['p_reference'] = 0.5 #bars! Just needs to be lower than your highest pressure

df = start_case.spectrum(opacity)

The other edit you will want to make is the addition of surface reflectance. For example, for a constant surface albedo of 0.3:

start_case.surface_reflect(0.3,opacity.wno) #this inputs the surface reflectivity and tells the code

This can be run anytime before you run your spectrum.

I will keep your question open until I put a proper raise Except. Thanks so much for pointing this out.

[Nicholaswogan]

Thanks! Looks like I had to add a slightly different line in v2.2.1

start_case.inputs['approx']['p_reference'] = 0.5

But ya everything looks good now!

[natashabatalha]

@Nicholaswogan I made a better change to just set reference pressure to the max pressure if it is out of bounds. For terrestrial planets especially, this makes more sense and results in the correct altitude calculations. Will merge to master and update pip to 2.3.1.

[Nicholaswogan]

Sounds great, thanks again.