Commit ace1dabb authored by Andrea Giannetti's avatar Andrea Giannetti
Browse files

Fixed output DB query multiplicity.

parent 0e620b9d
Loading
Loading
Loading
Loading
+1 −1
Original line number Original line Diff line number Diff line
@@ -410,7 +410,7 @@ def get_value_if_specified(parameters_dict: dict,
        return None
        return None




def parse_input_main() -> Tuple[str, str, float, float, list, int]:
def parse_input_main() -> Tuple[str, str, np.array, np.array, list, int]:
    """
    """
    Parse the stg and mdl configuration files, as well as the main one, to apply overrides, compute the arrays of
    Parse the stg and mdl configuration files, as well as the main one, to apply overrides, compute the arrays of
    temperatures and number densities and the type of models defined in the grid (isothermal/uniform density vs.
    temperatures and number densities and the type of models defined in the grid (isothermal/uniform density vs.
+1 −0
Original line number Original line Diff line number Diff line
@@ -210,6 +210,7 @@ def main_presentation_step(run_id: str,
            os.rmdir(scratches)
            os.rmdir(scratches)
    prs_inspection_main(run_id=run_id,
    prs_inspection_main(run_id=run_id,
                        is_isothermal=_tdust_model_type == 'isothermal',
                        is_isothermal=_tdust_model_type == 'isothermal',
                        is_homogeneous=_model_type == 'homogeneous',
                        engine=engine)
                        engine=engine)
    cleanup_tmp_table(run_id=run_id,
    cleanup_tmp_table(run_id=run_id,
                      engine=engine)
                      engine=engine)
+92 −45
Original line number Original line Diff line number Diff line
@@ -3,13 +3,14 @@ import numpy as np
import os
import os
import seaborn as sns
import seaborn as sns
import pandas as pd
import pandas as pd
import sqlalchemy
import xarray as xr
import xarray as xr
from typing import Union, Tuple, List
from typing import Union, Tuple, List
from contextlib import closing
from contextlib import closing
from itertools import product
from itertools import product
from astropy.io import fits
from astropy.io import fits
from astropy import units as u
from astropy import units as u
from sqlalchemy.orm import Session
from sqlalchemy.orm import Session, aliased
from sqlalchemy import and_, or_
from sqlalchemy import and_, or_
from sqlalchemy import engine as sqla_engine
from sqlalchemy import engine as sqla_engine
from stg.stg_build_db_structure import (GridPars,
from stg.stg_build_db_structure import (GridPars,
@@ -31,7 +32,9 @@ def get_aggregated_ratio_from_db(
        gas_density: float,
        gas_density: float,
        lines: Union[list, tuple],
        lines: Union[list, tuple],
        session: Session,
        session: Session,
        run_id: str) -> float:
        run_id: str,
        is_isothermal: bool = False,
        is_homogeneous: bool = False) -> float:
    """
    """
    Get the aggregated ratio from the database, according to the aggregation function specified in the pre config file;
    Get the aggregated ratio from the database, according to the aggregation function specified in the pre config file;
    this function works for a homogeneous model. For more complex models the query must be revised.
    this function works for a homogeneous model. For more complex models the query must be revised.
@@ -40,23 +43,50 @@ def get_aggregated_ratio_from_db(
    :param lines: the lines used to compute the ratio
    :param lines: the lines used to compute the ratio
    :param session: the SQLAlchemy session to use
    :param session: the SQLAlchemy session to use
    :param run_id: the id of the run
    :param run_id: the id of the run
    :param is_isothermal: whether the model is isothermal
    :param is_homogeneous: whether the model has a homogeneous density distribution
    :return: the aggregated ratio
    :return: the aggregated ratio
    """
    """
    results = session.query(GridPars, RatioMaps) \

        .join(StarsPars, isouter=True) \
    density_column, dust_temperature_column = get_filter_columns(is_homogeneous=is_homogeneous,
        .join(ModelPars) \
                                                                 is_isothermal=is_isothermal)
        .join(MomentZeroMaps) \
    mom_zero_1, mom_zero_2 = aliased(MomentZeroMaps), aliased(MomentZeroMaps)
        .join(RatioMaps,
    model_pars_1, model_pars_2 = aliased(ModelPars), aliased(ModelPars)
              or_(
    results = session.query(RatioMaps) \
                  and_(RatioMaps.mom_zero_map_1,
        .join(mom_zero_1,
                       ModelPars.iline == lines[0]),
              and_(RatioMaps.mom_zero_name_1 == mom_zero_1.mom_zero_name,
                  and_(RatioMaps.mom_zero_map_2,
                   RatioMaps.run_id == mom_zero_1.run_id)) \
                       ModelPars.iline == lines[1]))) \
        .join(model_pars_1,
              and_(mom_zero_1.fits_cube_name == model_pars_1.fits_cube_name,
                   mom_zero_1.run_id == model_pars_1.run_id)) \
        .join(mom_zero_2,
              and_(RatioMaps.mom_zero_name_2 == mom_zero_2.mom_zero_name,
                   RatioMaps.run_id == mom_zero_2.run_id)) \
        .join(model_pars_2,
              and_(mom_zero_2.fits_cube_name == model_pars_2.fits_cube_name,
                   mom_zero_2.run_id == model_pars_2.run_id)) \
        .join(GridPars) \
        .filter(
        .filter(
        and_(GridPars.dust_temperature_at_reference == dust_temperature,
        and_(dust_temperature_column == dust_temperature,
             GridPars.density_at_reference == gas_density,
             density_column == gas_density,
             GridPars.run_id == run_id)).order_by(GridPars.created_on.desc()).first()
             GridPars.run_id == run_id,
    return results[1].aggregated_ratio
             model_pars_1.iline == lines[0],
             model_pars_2.iline == lines[1])).order_by(GridPars.created_on.desc()).first()
    return results.aggregated_ratio


def get_filter_columns(is_homogeneous: bool,
                       is_isothermal: bool) -> Tuple[sqlalchemy.Column, sqlalchemy.Column]:
    # Select correct column for filtering
    if is_isothermal is True:
        dust_temperature_column = GridPars.dust_temperature
    else:
        dust_temperature_column = GridPars.dust_temperature_at_reference
    if is_homogeneous is True:
        density_column = GridPars.central_density
    else:
        density_column = GridPars.density_at_reference
    return density_column, dust_temperature_column




def get_results(
def get_results(
@@ -65,7 +95,8 @@ def get_results(
        lines: Union[list, tuple],
        lines: Union[list, tuple],
        run_id: str,
        run_id: str,
        session: Session,
        session: Session,
        is_isothermal: bool = False) -> Tuple[float, str, str, float, str, str, tuple, np.array]:
        is_isothermal: bool = False,
        is_homogeneous: bool = False) -> Tuple[float, str, str, float, str, str, tuple, np.array]:
    """
    """
    Get results from the database, given the parameters of the model
    Get results from the database, given the parameters of the model
    :param dust_temperature: the characteristic dust temperature
    :param dust_temperature: the characteristic dust temperature
@@ -74,33 +105,40 @@ def get_results(
    :param run_id: the run_id of the model
    :param run_id: the run_id of the model
    :param session: a session to query the database
    :param session: a session to query the database
    :param is_isothermal: whether the model is isothermal
    :param is_isothermal: whether the model is isothermal
    :param is_homogeneous: whether the model has a homogeneous density distribution
    :return: a tuple containing the aggregated line ratios, the moment zero maps names, the pixel size, the ratio map
    :return: a tuple containing the aggregated line ratios, the moment zero maps names, the pixel size, the ratio map
        name, the fits grid name (for the number density), the scaling factor between the grid and image pixel (must be
        name, the fits grid name (for the number density), the scaling factor between the grid and image pixel (must be
        integer!), and the pixel-by-pixel measured line ratios
        integer!), and the pixel-by-pixel measured line ratios
    """
    """
    # Select correct column for filtering
    density_column, dust_temperature_column = get_filter_columns(is_homogeneous=is_homogeneous,
    if is_isothermal is True:
                                                                 is_isothermal=is_isothermal)
        dust_temperature_column = GridPars.dust_temperature
    else:
        dust_temperature_column = GridPars.dust_temperature_at_reference


    mom_zero_1, mom_zero_2 = aliased(MomentZeroMaps), aliased(MomentZeroMaps)
    model_pars_1, model_pars_2 = aliased(ModelPars), aliased(ModelPars)
    # Query database
    # Query database
    results = session.query(GridPars, RatioMaps, GridFiles, ModelPars) \
    results = session.query(GridPars, RatioMaps, GridFiles, model_pars_1) \
        .join(mom_zero_1,
              and_(RatioMaps.mom_zero_name_1 == mom_zero_1.mom_zero_name,
                   RatioMaps.run_id == mom_zero_1.run_id)) \
        .join(model_pars_1,
              and_(mom_zero_1.fits_cube_name == model_pars_1.fits_cube_name,
                   mom_zero_1.run_id == model_pars_1.run_id)) \
        .join(mom_zero_2,
              and_(RatioMaps.mom_zero_name_2 == mom_zero_2.mom_zero_name,
                   RatioMaps.run_id == mom_zero_2.run_id)) \
        .join(model_pars_2,
              and_(mom_zero_2.fits_cube_name == model_pars_2.fits_cube_name,
                   mom_zero_2.run_id == model_pars_2.run_id)) \
        .join(GridPars) \
        .join(GridFiles) \
        .join(GridFiles) \
        .join(StarsPars, isouter=True) \
        .join(StarsPars, isouter=True) \
        .join(ModelPars) \
        .join(MomentZeroMaps) \
        .join(RatioMaps,
              or_(
                  and_(RatioMaps.mom_zero_map_1,
                       ModelPars.iline == lines[0]),
                  and_(RatioMaps.mom_zero_map_2,
                       ModelPars.iline == lines[1]))) \
        .filter(
        .filter(
        and_(dust_temperature_column == dust_temperature,
        and_(dust_temperature_column == dust_temperature,
             GridPars.density_at_reference == gas_density,
             density_column == gas_density,
             GridFiles.quantity == 'gas_number_density',
             GridFiles.quantity == 'gas_number_density',
             GridPars.run_id == run_id)).order_by(GridPars.created_on.desc()).first()
             GridPars.run_id == run_id),
        model_pars_1.iline == lines[0],
        model_pars_2.iline == lines[1]).order_by(GridPars.created_on.desc()).first()
    with fits.open(os.path.join('prs', 'fits', 'ratios', results[1].ratio_map_name)) as ratio_fitsfile:
    with fits.open(os.path.join('prs', 'fits', 'ratios', results[1].ratio_map_name)) as ratio_fitsfile:
        ratio_values = ratio_fitsfile[0].data
        ratio_values = ratio_fitsfile[0].data
    assert (results[3].npix % results[0].grid_shape_1 == 0) and (results[3].npix % results[0].grid_shape_2 == 0)
    assert (results[3].npix % results[0].grid_shape_1 == 0) and (results[3].npix % results[0].grid_shape_2 == 0)
@@ -120,7 +158,8 @@ def get_grid_values(
        gas_density: float,
        gas_density: float,
        run_id: str,
        run_id: str,
        session: Session,
        session: Session,
        is_isothermal: bool = False) -> np.array:
        is_isothermal: bool = False,
        is_homogeneous: bool = False) -> np.array:
    """
    """
    Retrieve the grid values, given the model parameters
    Retrieve the grid values, given the model parameters
    :param quantity_name: the quantity to extract
    :param quantity_name: the quantity to extract
@@ -129,20 +168,18 @@ def get_grid_values(
    :param run_id: the run_id to process
    :param run_id: the run_id to process
    :param session: a session to query the database
    :param session: a session to query the database
    :param is_isothermal: whether the model is isothermal
    :param is_isothermal: whether the model is isothermal
    :param is_homogeneous: whether the model has a homogeneous density distribution
    :return: the array of grid values
    :return: the array of grid values
    """
    """
    # Select correct column for filtering
    density_column, dust_temperature_column = get_filter_columns(is_homogeneous=is_homogeneous,
    if is_isothermal is True:
                                                                 is_isothermal=is_isothermal)
        dust_temperature_column = GridPars.dust_temperature
    else:
        dust_temperature_column = GridPars.dust_temperature_at_reference


    # Query database
    # Query database
    results = session.query(GridPars, GridFiles) \
    results = session.query(GridPars, GridFiles) \
        .join(GridFiles) \
        .join(GridFiles) \
        .filter(
        .filter(
        and_(dust_temperature_column == dust_temperature,
        and_(dust_temperature_column == dust_temperature,
             GridPars.density_at_reference == gas_density,
             density_column == gas_density,
             GridFiles.quantity == quantity_name,
             GridFiles.quantity == quantity_name,
             GridPars.run_id == run_id)).order_by(GridPars.created_on.desc()).first()
             GridPars.run_id == run_id)).order_by(GridPars.created_on.desc()).first()
    with fits.open(os.path.join('prs', 'fits', 'grids', results[1].fits_grid_name)) as fitsfile:
    with fits.open(os.path.join('prs', 'fits', 'grids', results[1].fits_grid_name)) as fitsfile:
@@ -171,6 +208,7 @@ def get_column_density_vs_mom0(volume_density_grid: np.array,


def main(run_id: str,
def main(run_id: str,
         is_isothermal: bool,
         is_isothermal: bool,
         is_homogeneous: bool,
         engine: sqla_engine = None):
         engine: sqla_engine = None):
    if engine is None:
    if engine is None:
        engine = get_pg_engine(logger=logger)
        engine = get_pg_engine(logger=logger)
@@ -199,7 +237,9 @@ def main(run_id: str,
                                                                gas_density=nh2,
                                                                gas_density=nh2,
                                                                lines=lines,
                                                                lines=lines,
                                                                session=session,
                                                                session=session,
                                                                run_id=run_id)
                                                                run_id=run_id,
                                                                is_isothermal=is_isothermal,
                                                                is_homogeneous=is_homogeneous)
                results.loc[tdust, nh2] = aggregated_ratio
                results.loc[tdust, nh2] = aggregated_ratio
                logger.debug(f'The aggregated ratio for lines {lines}, using {nh2}, {tdust} is: {aggregated_ratio}')
                logger.debug(f'The aggregated ratio for lines {lines}, using {nh2}, {tdust} is: {aggregated_ratio}')
                ratio, \
                ratio, \
@@ -216,7 +256,8 @@ def main(run_id: str,
                        lines=lines,
                        lines=lines,
                        session=session,
                        session=session,
                        run_id=run_id,
                        run_id=run_id,
                        is_isothermal=is_isothermal
                        is_isothermal=is_isothermal,
                        is_homogeneous=is_homogeneous
                    )
                    )
                density_grid = get_grid_values(
                density_grid = get_grid_values(
                    quantity_name='gas_number_density',
                    quantity_name='gas_number_density',
@@ -224,7 +265,8 @@ def main(run_id: str,
                    gas_density=nh2,
                    gas_density=nh2,
                    session=session,
                    session=session,
                    run_id=run_id,
                    run_id=run_id,
                    is_isothermal=is_isothermal
                    is_isothermal=is_isothermal,
                    is_homogeneous=is_homogeneous
                )
                )
                temperature_grid = get_grid_values(
                temperature_grid = get_grid_values(
                    quantity_name='dust_temperature',
                    quantity_name='dust_temperature',
@@ -232,7 +274,8 @@ def main(run_id: str,
                    gas_density=nh2,
                    gas_density=nh2,
                    session=session,
                    session=session,
                    run_id=run_id,
                    run_id=run_id,
                    is_isothermal=is_isothermal
                    is_isothermal=is_isothermal,
                    is_homogeneous=is_homogeneous
                )
                )
                avg_density_map = np.nansum(np.where(density_grid == 0, np.nan, density_grid ** 2), axis=2) / \
                avg_density_map = np.nansum(np.where(density_grid == 0, np.nan, density_grid ** 2), axis=2) / \
                                  np.nansum(np.where(density_grid == 0, np.nan, density_grid), axis=2)
                                  np.nansum(np.where(density_grid == 0, np.nan, density_grid), axis=2)
@@ -279,6 +322,10 @@ def main(run_id: str,
        sns.kdeplot(x=np.log10(df_all['avg_nh2']), y=df_all[f'ratio_{"-".join(lines)}'], bw_adjust=2)
        sns.kdeplot(x=np.log10(df_all['avg_nh2']), y=df_all[f'ratio_{"-".join(lines)}'], bw_adjust=2)
        plt.savefig(os.path.join('prs', 'output', f'ratio_vs_avg_density_los_kde_{"-".join(lines)}.png'))
        plt.savefig(os.path.join('prs', 'output', f'ratio_vs_avg_density_los_kde_{"-".join(lines)}.png'))
        plt.clf()
        plt.clf()
    # for lines in line_pairs:
    #     sns.kdeplot(x=np.log10(df_all['std_nh2']), y=df_all[f'ratio_{"-".join(lines)}'], bw_adjust=2)
    #     plt.savefig(os.path.join('prs', 'output', f'ratio_vs_std_density_los_kde_{"-".join(lines)}.png'))
    #     plt.clf()
    df_all.to_csv('full_dataset.csv')
    df_all.to_csv('full_dataset.csv')




@@ -311,5 +358,5 @@ def plot_ratio_vs_density(lines: Tuple[str, str],




if __name__ == '__main__':
if __name__ == '__main__':
    main(run_id='0d0761f8-fcf3-4293-9d9d-ba1a37f79a19', is_isothermal=False)
    main(run_id='0d0761f8-fcf3-4293-9d9d-ba1a37f79a19', is_isothermal=False, is_homogeneous=False)
    # main(run_id='ba7fd3ed-5947-4dc6-bcef-38f151f19b77', is_isothermal=False)
    # main(run_id='ba7fd3ed-5947-4dc6-bcef-38f151f19b77', is_isothermal=False)