Telescope class and generalization

.gitignore

@@ -1,5 +1,6 @@
 *.pyc
 *.iml
+*.sh
 .ipynb_checkpoints/
 config_local.ini
 .idea/

pastis/config.py

@@ -20,9 +20,12 @@ def get_config_ini_path():
     return os.path.join(os.path.dirname(os.path.realpath(__file__)), config_file_name)
 
 
-def load_config_ini():
+def load_config_ini(dir=None):
     # Locate where on disk this file is.
-    code_directory = os.path.dirname(os.path.realpath(__file__))
+    if dir == None:
+        code_directory = os.path.dirname(os.path.realpath(__file__))
+    else:
+        code_directory = dir
 
     # Read config file once here.
     config = configparser.ConfigParser()

pastis/config_pastis.ini

@@ -6,6 +6,9 @@
 ; figure out webbpsf-data path with: webbpsf.utils.get_webbpsf_data_path()
 webbpsf_data_path = /Users/<user-name>/anaconda/envs/astroconda/share/webbpsf-data
 local_data_path = /Users/<user-name>/data_from_repos/pastis_data
+analysis_name = last
+; folder name to analysis or hockeystick curve, or 'last' to run last folder.
+
 
 [telescope]
 name = LUVOIR
@@ -131,6 +134,7 @@ valid_range_lower = -4
 valid_range_upper = 4
 
 ; telescope
+design = small
 nb_subapertures = 120
 diameter = 15.
 gaps = 0.02
@@ -141,6 +145,16 @@ lyot_stop_path_in_optics = inputs/LS_LUVOIR_ID0120_OD0982_no_struts_gy_ovsamp4_N
 
 ; absolute path to Harris spreadsheet, currently hosted off the repo
 harris_data_path = /Users/yourname/somewhere/Sensitivities2.xlsx
+therm = False
+mech = True
+other = False
+
+; Possible: "seg_mirror", "harris_seg_mirror", "zernike_mirror", "default"
+DM = zernike_mirror
+; 1 will be set for harris_seg_mirror, "default" limit at 3
+DM_mode = 1
+; 0 to do all mode, accept multiple values but all must be lower than DM_mode
+DM_mode_select = 0
 
 ; coronagraph
 ; iwa and owa from dictionaries within files. could move that to util.
@@ -191,6 +205,31 @@ IWA = 2
 OWA = 32
 lyot_stop_ratio = 0.95
 
+[generation]
+run = True
+
+method = intensity
+;Valid method: intensity, Efield
+
+[hockeystick]
+run = True
+
+range_points = 30
+no_realizations = 1
+
+[analysis]
+run = True
+
+calculate_modes = True
+calculate_sigmas = True
+run_monte_carlo_modes = True
+calc_cumulative_contrast = True
+calculate_mus = True
+run_monte_carlo_segments = True
+calculate_covariance_matrices = True
+analytical_statistics = True
+calculate_segment_based = True
+
 [numerical]
 ; size_seg used to be 100 in atlast case, 118 for JWST 512 px images, 239 for JWST 1024 px images
 size_seg = 239
@@ -206,6 +245,7 @@ z_pup_downsample = 10
 
 ; this is not used automatically in the functions, it is always defined (or read from here) manually
 current_analysis = 2020-01-13T21-34-29_luvoir-small
+; deprecated use [local] analysis_name
 
 [zernikes]
 ; Noll convention!
@@ -234,3 +274,10 @@ number_of_low_order_modes = 15
 number_of_mid_order_modes = 1
 number_of_high_order_modes = 4
 number_of_continuous_dm_actuators = 4
+
+[save_data]
+save_efields = True
+save_psfs = True
+save_opds = True
+save_coro_floor = True
+coro_simulator = True

pastis/contrast_calculation_simple.py

@@ -20,6 +20,7 @@
 from pastis.config import CONFIG_PASTIS
 from pastis.e2e_simulators.hicat_imaging import set_up_hicat
 from pastis.e2e_simulators.luvoir_imaging import LuvoirAPLC
+from pastis.launchers.parameters import parameters
 import pastis.e2e_simulators.webbpsf_imaging as webbpsf_imaging
 import pastis.analytical_pastis.image_pastis as impastis
 import pastis.util as util
@@ -414,15 +415,14 @@ def contrast_rst_num(coro_floor, norm, matrix_dir, rms=50*u.nm):
     # Create random aberration coefficients on segments, scaled to total rms
     aber = util.create_random_rms_values(total_seg, rms)
 
-    ### E2E JWST sim
+    ### E2E RST sim
     start_e2e = time.time()
 
     rst_sim = webbpsf_imaging.set_up_cgi()
     rst_sim.fpm = CONFIG_PASTIS.get('RST', 'fpm')
     nb_actu = rst_sim.nbactuator
     iwa = CONFIG_PASTIS.getfloat('RST', 'IWA')
     owa = CONFIG_PASTIS.getfloat('RST', 'OWA')
-    sampling = CONFIG_PASTIS.getfloat('RST', 'sampling')
 
     # Put aberration on OTE
     rst_sim.dm1.flatten()
@@ -462,6 +462,67 @@ def contrast_rst_num(coro_floor, norm, matrix_dir, rms=50*u.nm):
     return contrast_rst, contrast_matrix
 
 
+def contrast_general_num(matrix_dir, rms=50*u.nm):
+    """
+    Compute the contrast for a random aberration over all DM actuators in the RST simulator.
+
+    :param matrix_dir: str, directory of saved matrix
+    :param rms: astropy quantity (e.g. m or nm), WFE rms (OPD) to be put randomly over the entire continuous mirror
+    :return: 2x float, E2E and matrix contrast
+    """
+    # Keep track of time
+    start_time = time.time()
+
+    # Parameters
+    param = parameters()
+    telescope = param.def_telescope()
+    telescope.normalization_and_dark_hole()
+    telescope.calculate_unaberrated_contrast()
+    total_seg = telescope.number_all_modes
+
+    # Import numerical PASTIS matrix
+    filename = 'pastis_matrix'
+    matrix_pastis = fits.getdata(os.path.join(matrix_dir, filename + '.fits'))
+
+    # Create random aberration coefficients on segments, scaled to total rms
+    aber = util.create_random_rms_values(total_seg, rms)
+
+    ### E2E RST sim
+    start_e2e = time.time()
+
+    # Put aberration on OTE
+    telescope.flatten()
+    for nseg in range(total_seg):
+        telescope.push_mode(nseg, aber[nseg].value*u.nm)
+
+    telescope.imaging_psf()
+
+    # Get the mean contrast
+    contrast = telescope.contrast()
+    end_e2e = time.time()
+
+    ## MATRIX PASTIS
+    log.info('Generating contrast from matrix-PASTIS')
+    start_matrixpastis = time.time()
+    # Get mean contrast from matrix PASTIS
+    contrast_matrix = util.pastis_contrast(aber, matrix_pastis) + telescope.contrast_floor   # calculating contrast with PASTIS matrix model
+    end_matrixpastis = time.time()
+
+    ## Outputs
+    log.info('\n--- CONTRASTS: ---')
+    log.info(f'Mean contrast from E2E: {contrast}')
+    log.info(f'Contrast from matrix PASTIS: {contrast_matrix}')
+
+    log.info('\n--- RUNTIMES: ---')
+    log.info(f'E2E: {end_e2e-start_e2e}sec = {(end_e2e-start_e2e)/60}min')
+    log.info(f'Matrix PASTIS: {end_matrixpastis-start_matrixpastis}sec = {(end_matrixpastis-start_matrixpastis)/60}min')
+
+    end_time = time.time()
+    runtime = end_time - start_time
+    log.info(f'Runtime for contrast_calculation_simple.py: {runtime} sec = {runtime/60} min')
+
+    return contrast, contrast_matrix
+
 if __name__ == '__main__':
 
     # Test JWST