<commit_before>
import os
from spectral_cube import SpectralCube
from spectral_cube.lower_dimensional_structures import Projection
import astropy.units as u
from astropy.coordinates import SkyCoord
import numpy as np
from astropy.io import fits
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as p

from basics.bubble_segment3D import BubbleFinder
from basics.utils import sig_clip

'''
Create the bubble catalogue of M33 with the 14B-088 map
'''

data_path = "/media/eric/MyRAID/M33/14B-088/HI/full_imaging/"
# data_path = "/lustre/home/ekoch/m33/"

cube = \
    SpectralCube.read(os.path.join(data_path,
                                    "M33_14B-088_HI.clean.image.pbcov_gt_0.3_masked.fits"))

mom2 = fits.getdata(os.path.join(data_path,
    "M33_14B-088_HI.clean.image.pbcov_gt_0.3_masked.rotsub.lwidth.fits"))

sigma = sig_clip(cube[-1].value, nsig=10)

# Create a cruddy linewidth map. This will need to be improved, but that can
# be done with just the saved bubble objects
# lwidth = cube.with_mask(cube > 2 * sigma * u.Jy).linewidth_sigma()
lwidth = np.sqrt(mom2) * u.m / u.s

galaxy_props = {"center_coord":
                SkyCoord(23.461667, 30.660194,
                         unit=(u.deg, u.deg), frame='fk5'),
                "inclination": 56. * u.deg,
                "position_angle": 201. * u.deg,
                "scale_height": 100. * u.pc}


scales = 3. * np.arange(1, 15, np.sqrt(2))

# The first few channels have some emission in them. It's partially M33,
# partially galactic HI. So estimate the noise level from the last channel
# which is pretty much noise only. It gives 1.8 mJy/bm, which is right on.
bub_find = BubbleFinder(cube, keep_threshold_mask=True,
                        # empty_channel=cube.shape[0] - 1,  # Given above
                        sigma=sigma,
                        galaxy_props=galaxy_props,
                        distance=0.84 * u.Mpc)

# Setup folder for saving the 2D regions
folder = os.path.join(data_path, "bubbles_2D")
try:
    os.mkdir(folder)
except OSError:
    pass

# Requiring minimum of 15 channels, similar to the 3 channel requirement with
# THINGS (15 * 0.2 km/s = 3 km/s; 3 * 1.2 km/s = 3.6 km/s)
bub_find.get_bubbles(verbose=True, overlap_frac=0.5, multiprocess=True,
                     refit=False, nsig=1.5, min_corr=0.7, min_overlap=0.7,
                     global_corr=0.5, min_channels=15, nprocesses=None,
                     scales=scales, cube_linewidth=lwidth, save_regions=True,
                     save_region_path=os.path.join(data_path, "bubbles_2D"))

# Save all of the bubbles
# folder = os.path.expanduser("~/MyRAID/M33/14B-088/HI/full_imaging/bubbles/")
folder = os.path.join(data_path, "bubbles/")
try:
    os.mkdir(folder)
except OSError:
    pass

name = "M33_14B-088"

bub_find.save_bubbles(folder=folder, name=name)

catalog = bub_find.to_catalog()

catalog.write_table(os.path.join(folder, "bubble_catalog.ecsv"))

# Save the mask
np.save(os.path.join(folder, "{}_bubble_mask.npy".format(name)),
        bub_find.mask)


# Show the outlines of the bubbles
fig = p.figure()
ax = fig.add_subplot(111)
ax = bub_find.visualize_bubbles(show=False, ax=ax,
                                plot_twoD_shapes=False)
fig.savefig(os.path.join(folder, "{}_mom0_bubbles.pdf".format(name)))
p.close()
# With the 2D regions
fig = p.figure()
ax = fig.add_subplot(111)
ax = bub_find.visualize_bubbles(show=False, ax=ax,
                                plot_twoD_shapes=True)
fig.savefig(os.path.join(folder,
                         "{}_mom0_bubbles_w_twoD.pdf".format(name)))
p.close()

# Individual channel maps
output_folder_chans = os.path.join(folder, "channel_maps")
try:
    os.mkdir(output_folder_chans)
except OSError:
    pass
bub_find.visualize_channel_maps(show_mask_contours=True, plot_unclustered=True,
                                save=True, save_path=output_folder_chans,
                                save_name=name)
<commit_msg>Loading in the linewidth, not moment 2...<commit_after>
import os
from spectral_cube import SpectralCube
from spectral_cube.lower_dimensional_structures import Projection
import astropy.units as u
from astropy.coordinates import SkyCoord
import numpy as np
from astropy.io import fits
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as p

from basics.bubble_segment3D import BubbleFinder
from basics.utils import sig_clip

'''
Create the bubble catalogue of M33 with the 14B-088 map
'''

data_path = "/media/eric/MyRAID/M33/14B-088/HI/full_imaging/"
# data_path = "/lustre/home/ekoch/m33/"

cube = \
    SpectralCube.read(os.path.join(data_path,
                                    "M33_14B-088_HI.clean.image.pbcov_gt_0.3_masked.fits"))

lwidth = fits.getdata(os.path.join(data_path,
    "M33_14B-088_HI.clean.image.pbcov_gt_0.3_masked.rotsub.lwidth.fits"))

sigma = sig_clip(cube[-1].value, nsig=10)

# Create a cruddy linewidth map. This will need to be improved, but that can
# be done with just the saved bubble objects
# lwidth = cube.with_mask(cube > 2 * sigma * u.Jy).linewidth_sigma()
lwidth = lwidth * u.m / u.s

galaxy_props = {"center_coord":
                SkyCoord(23.461667, 30.660194,
                         unit=(u.deg, u.deg), frame='fk5'),
                "inclination": 56. * u.deg,
                "position_angle": 201. * u.deg,
                "scale_height": 100. * u.pc}


scales = 3. * np.arange(1, 15, np.sqrt(2))

# The first few channels have some emission in them. It's partially M33,
# partially galactic HI. So estimate the noise level from the last channel
# which is pretty much noise only. It gives 1.8 mJy/bm, which is right on.
bub_find = BubbleFinder(cube, keep_threshold_mask=True,
                        # empty_channel=cube.shape[0] - 1,  # Given above
                        sigma=sigma,
                        galaxy_props=galaxy_props,
                        distance=0.84 * u.Mpc)

# Setup folder for saving the 2D regions
folder = os.path.join(data_path, "bubbles_2D")
try:
    os.mkdir(folder)
except OSError:
    pass

# Requiring minimum of 15 channels, similar to the 3 channel requirement with
# THINGS (15 * 0.2 km/s = 3 km/s; 3 * 1.2 km/s = 3.6 km/s)
bub_find.get_bubbles(verbose=True, overlap_frac=0.5, multiprocess=True,
                     refit=False, nsig=1.5, min_corr=0.7, min_overlap=0.7,
                     global_corr=0.5, min_channels=15, nprocesses=None,
                     scales=scales, cube_linewidth=lwidth, save_regions=True,
                     save_region_path=os.path.join(data_path, "bubbles_2D"))

# Save all of the bubbles
# folder = os.path.expanduser("~/MyRAID/M33/14B-088/HI/full_imaging/bubbles/")
folder = os.path.join(data_path, "bubbles/")
try:
    os.mkdir(folder)
except OSError:
    pass

name = "M33_14B-088"

bub_find.save_bubbles(folder=folder, name=name)

catalog = bub_find.to_catalog()

catalog.write_table(os.path.join(folder, "bubble_catalog.ecsv"))

# Save the mask
np.save(os.path.join(folder, "{}_bubble_mask.npy".format(name)),
        bub_find.mask)


# Show the outlines of the bubbles
fig = p.figure()
ax = fig.add_subplot(111)
ax = bub_find.visualize_bubbles(show=False, ax=ax,
                                plot_twoD_shapes=False)
fig.savefig(os.path.join(folder, "{}_mom0_bubbles.pdf".format(name)))
p.close()
# With the 2D regions
fig = p.figure()
ax = fig.add_subplot(111)
ax = bub_find.visualize_bubbles(show=False, ax=ax,
                                plot_twoD_shapes=True)
fig.savefig(os.path.join(folder,
                         "{}_mom0_bubbles_w_twoD.pdf".format(name)))
p.close()

# Individual channel maps
output_folder_chans = os.path.join(folder, "channel_maps")
try:
    os.mkdir(output_folder_chans)
except OSError:
    pass
bub_find.visualize_channel_maps(show_mask_contours=True, plot_unclustered=True,
                                save=True, save_path=output_folder_chans,
                                save_name=name)

{'fraction_non_alphanumeric': 0.11349036402569593, 'fraction_numerical': 0.032875677037410254, 'mean_word_length': 2.8788470933072787, 'pattern_counts': {'":': 8, '<': 3, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '344911', 'n_tokens_mistral': 3008, 'n_tokens_neox': 2762, 'n_words': 598}