From 826e74efcecc666fa3848fcec62009607e88d615 Mon Sep 17 00:00:00 2001 From: Carlo Castoldi Date: Mon, 25 Nov 2024 16:59:37 +0100 Subject: [PATCH 1/5] initial work to add Kim DevCCFv1 2024 --- .../atlas_scripts/kim_devccf_mouse.py | 258 ++++++++++++++++++ 1 file changed, 258 insertions(+) create mode 100644 brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py diff --git a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py new file mode 100644 index 00000000..baac91db --- /dev/null +++ b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py @@ -0,0 +1,258 @@ +import json +import multiprocessing as mp +import time +from pathlib import Path + +import numpy as np +import pandas as pd +from brainglobe_utils.IO.image import load_nii +from rich.progress import track + +from brainglobe_atlasapi.atlas_generation.mesh_utils import ( + Region, + create_region_mesh, +) +from brainglobe_atlasapi.atlas_generation.wrapup import wrapup_atlas_from_data + +# from brainglobe_atlasapi.config import DEFAULT_WORKDIR +from brainglobe_atlasapi.structure_tree_util import get_structures_tree + +ATLAS_NAME = "kim_devccf_mouse" # multiple versions of the same atlas +SPECIES = "Mus musculus" +AGE = "P14" +ATLAS_LINK = "https://kimlab.io/brain-map/DevCCF/" +CITATION = "Kronman, F.N., Liwang, J.K., Betty, R. et al. 2024, https://doi.org/10.1038/s41467-024-53254-w" +ORIENTATION = "posterior superior left" +ROOT_ID = 15564 +RESOLUTION_UM = 20 +VERSION = 1 +PACKAGER = "Carlo Castoldi " +ATLAS_FILE_URL = "https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/download.aspx?UniqueId=fe3d1692%2D94e4%2D4238%2Db6bc%2D95d18bcac022" +# curl 'https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/download.aspx?UniqueId=fe3d1692%2D94e4%2D4238%2Db6bc%2D95d18bcac022' -H 'User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:132.0) Gecko/20100101 Firefox/132.0' -H 'Accept: text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8' -H 'Accept-Language: en-GB,en;q=0.5' -H 'Accept-Encoding: gzip, deflate, br, zstd' -H 'Referer: https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/onedrive.aspx?id=%2Fpersonal%2Fyuk17%5Fpsu%5Fedu%2FDocuments%2FProjects%2FDevelopmental%5FAtlasing%2FDevCCF%5FTeamFileshare%2FDevCCF%5FTemplates%2FDevCCFv1%2B%2FDevCCFv1' -H 'Upgrade-Insecure-Requests: 1' -H 'Sec-Fetch-Dest: iframe' -H 'Sec-Fetch-Mode: navigate' -H 'Sec-Fetch-Site: same-origin' -H 'Connection: keep-alive' -H 'Cookie: FedAuth=77u/PD94bWwgdmVyc2lvbj0iMS4wIiBlbmNvZGluZz0idXRmLTgiPz48U1A+VjEzLDBoLmZ8bWVtYmVyc2hpcHx1cm4lM2FzcG8lM2Fhbm9uIzc1Yjc1NjA0NGUxZDI5ZDU5NjFmMjZjODc0MzFlZmQwNmU4ZWM2NzYwMTJmZGEyNjgxYzg4Zjc5ODdmMGY2YTAsMCMuZnxtZW1iZXJzaGlwfHVybiUzYXNwbyUzYWFub24jNzViNzU2MDQ0ZTFkMjlkNTk2MWYyNmM4NzQzMWVmZDA2ZThlYzY3NjAxMmZkYTI2ODFjODhmNzk4N2YwZjZhMCwxMzM3NjQ5NTA5NzAwMDAwMDAsMCwxMzM3NjY3MzA4NzYxOTQyMTMsMC4wLjAuMCwyNTgsN2NmNDhkNDUtM2RkYi00Mzg5LWE5YzEtYzExNTUyNmViNTJlLCwsZWU2YTY1YTEtMDBmNS02MDAwLWZiMDYtYTgzYjk1Yjg3NTIzLDkxYzI2NWExLTMwNDctNjAwMC1mYjA2LWFiY2IwMTQ2MTZmNSxDR1EvcVlGcGJFeVlJTHAvM2NHb3pnLDAsMCwwLCwsLDI2NTA0Njc3NDM5OTk5OTk5OTksMCwsLCwsLCwwLCwxOTI2MDksdVhlaFFKUGxlVmpOQ2Jha1VoR0Q2SXlGUVFrLFkzdXYyY1dWdUZMWElneHRSOW4vMVVValJ1UUNMd2wybzdZbE51S0RtNTBxOE16ODJVaCt6K0ZTL0NqUU0rNlZrc090RlkwY211OUJNWHllc2U3Z3dST1VMcXRDWEorV3M1NkVzcEZrdk1qTEI0VzNhbDVZeVRlZ0greURRaUFRbG8rVG82aUVmVndVODU5Z2RWME5GSUprenJSVUlCeWVDUThiMmhMRjhQNmlJTDZHYldBNXVwRnJpcTl6cml6THhheFBRUWtVRzRoUmJjMEgzaVMzMnVKekw1ckFuWnVZSFk1N2J4K2xZOGE5MUNVQ3BSQXJ3aU45NHZ5aXdaQWdlcDRXUkxwWlVWOU5NNFZwZ0xUN05zVUh2K0JBcGwxbnViY0FDVENlWElhbWNwclZ0eTBFdElpajVSeXR3UGk0OFRlSVlIQXVrRjhDUFRFVEdDcTlGQT09PC9TUD4=; FeatureOverrides_experiments=[]; MicrosoftApplicationsTelemetryDeviceId=8293cff3-d815-447c-b343-546e678a0450; MSFPC=GUID=79c2671112584c09b960ec538dca634d&HASH=79c2&LV=202411&V=4&LU=1732022021463; ai_session=O/d7EbWPiaNaTZQh0qUYpP|1732113248412|1732113248412' --output P14.zip +PARALLEL = False # disable parallel mesh extraction for easier debugging + + +TIMEPOINTS = ( + "E11.5", + "E13.5", + "E15.5", + "E18.5", + "P04", + "P14", + "P56" +) +RESOLUTIONS = (20, 50) +MODALITIES = ( + "LSFM", # Light Sheet Fluorescence Microscopy + "MRI-adc", # MRI Apparent Diffusion Coefficient + "MRI-dwi", # MRI Difusion Weighted Imaging + "MRI-fa", # MRI Fractional Anisotropy + "MRI-MTR", # MRI Magnetization Transfer Ratio + "MRI-T2" # MRI T2-weighted +) + +def get_annotations(devccfv1_root: str|Path, + age: str, + resolution: int): + if not isinstance(devccfv1_root, Path): + devccfv1_root = Path(devccfv1_root) + assert age in TIMEPOINTS, f"Unknown age timepoint: '{age}'" + assert resolution in RESOLUTIONS, f"Unknown resolution in µm: '{resolution}'" + annotations_path = devccfv1_root/age/f"{age}_DevCCF_Annotations_{resolution}um.nii.gz" + annotations = load_nii(annotations_path, as_array=True) + return annotations + +def get_reference(devccfv1_root: str|Path, + age: str, + resolution: int, + modality: str): + if not isinstance(devccfv1_root, Path): + devccfv1_root = Path(devccfv1_root) + assert age in TIMEPOINTS, f"Unknown age timepoint: '{age}'" + assert resolution in RESOLUTIONS, f"Unknown resolution: '{resolution}'" + assert modality in MODALITIES, f"Unknown modality: '{modality}'" + reference_path = devccfv1_root/age/f"{age}_{modality}_{resolution}um.nii.gz" + reference = load_nii(reference_path, as_array=True) + return reference + +def get_ontology(devccfv1_root: str|Path): + if not isinstance(devccfv1_root, Path): + devccfv1_root = Path(devccfv1_root) + DevCCFv1_path = devccfv1_root/"DevCCFv1_OntologyStructure.xlsx" + xl = pd.ExcelFile(DevCCFv1_path) + # xl.sheet_names # it has two excel sheets + # 'DevCCFv1_Ontology', 'README' + df = xl.parse("DevCCFv1_Ontology", header=1) + df = df[["Acronym", "ID16", "Name", "Structure ID Path16", "R", "G", "B"]] + df.rename(columns={ + "Acronym": "acronym", + "ID16": "id", + "Name": "name", + "Structure ID Path16": "structure_id_path", + "R": "r", + "G": "g", + "B": "b" + }, inplace=True) + structures = list(df.to_dict(orient="index").values()) + for structure in structures: + if structure["acronym"] == "mouse": + structure["acronym"] = "root" + structure_path = structure["structure_id_path"] + structure["structure_id_path"] = [int(id) for id in structure_path.strip("/").split("/")] + structure["rgb_triplet"] = [structure["r"], structure["g"], structure["b"]] + del structure["r"] + del structure["g"] + del structure["b"] + return structures + +def create_meshes(download_dir_path: str|Path, + structures, annotation_volume, root_id): + if not isinstance(download_dir_path, Path): + download_dir_path = Path(download_dir_path) + meshes_dir_path = download_dir_path / "meshes" + meshes_dir_path.mkdir(exist_ok=True) + + tree = get_structures_tree(structures) + + labels = np.unique(annotation_volume).astype(np.uint16) + + for key, node in tree.nodes.items(): + if key in labels: + is_label = True + else: + is_label = False + node.data = Region(is_label) + + # Mesh creation + closing_n_iters = 2 + decimate_fraction = 0.2 # 0.04 + # What fraction of the original number of vertices is to be kept. + smooth = False # smooth meshes after creation + start = time.time() + if PARALLEL: + pool = mp.Pool(mp.cpu_count() - 2) + try: + pool.map( + create_region_mesh, + [ + ( + meshes_dir_path, + node, + tree, + labels, + annotation_volume, + root_id, + closing_n_iters, + decimate_fraction, + smooth, + ) + for node in tree.nodes.values() + ], + ) + except mp.pool.MaybeEncodingError: + pass + else: + for node in track( + tree.nodes.values(), + total=tree.size(), + description="Creating meshes", + ): + # root_node = tree.nodes[root_id] + create_region_mesh( + ( + meshes_dir_path, + # root_node, + node, + tree, + labels, + annotation_volume, + root_id, + closing_n_iters, + decimate_fraction, + smooth, + ) + ) + + print( + "Finished mesh extraction in: ", + round((time.time() - start) / 60, 2), + " minutes", + ) + return meshes_dir_path + + +def create_mesh_dict(structures, meshes_dir_path): + meshes_dict = dict() + structures_with_mesh = [] + for s in structures: + # Check if a mesh was created + mesh_path = meshes_dir_path / f'{s["id"]}.obj' + if not mesh_path.exists(): + print(f"No mesh file exists for: {s}, ignoring it") + continue + else: + # Check that the mesh actually exists (i.e. not empty) + if mesh_path.stat().st_size < 512: + print(f"obj file for {s} is too small, ignoring it.") + continue + + structures_with_mesh.append(s) + meshes_dict[s["id"]] = mesh_path + + print( + f"In the end, {len(structures_with_mesh)} " + "structures with mesh are kept" + ) + return meshes_dict, structures_with_mesh + + +if __name__ == "__main__": + atlas_root = "/home/castoldi/Downloads/DevCCFv1" + # Generated atlas path: + DEFAULT_WORKDIR = Path.home()/"brainglobe_workingdir" + bg_root_dir = DEFAULT_WORKDIR/ATLAS_NAME + download_dir_path = bg_root_dir/"downloads" + download_dir_path.mkdir(exist_ok=True, parents=True) + + structures = get_ontology(atlas_root) + # save regions list json: + with open(download_dir_path/"structures.json", "w") as f: + json.dump(structures, f) + annotation_volume = get_annotations(atlas_root, AGE, RESOLUTION_UM) + reference_volume = get_reference(atlas_root, AGE, RESOLUTION_UM, "LSFM") + # Create meshes: + print(f"Saving atlas data at {download_dir_path}") + meshes_dir_path = create_meshes( + download_dir_path, + structures, + annotation_volume, + ROOT_ID + ) + + meshes_dict, structures_with_mesh = create_mesh_dict( + structures, meshes_dir_path + ) + # Wrap up, compress, and remove file: + print("Finalising atlas") + output_filename = wrapup_atlas_from_data( + atlas_name=ATLAS_NAME, + atlas_minor_version=VERSION, + citation=CITATION, + atlas_link=ATLAS_LINK, + species=SPECIES, + resolution=RESOLUTION_UM, + orientation=ORIENTATION, + root_id=ROOT_ID, + reference_stack=reference_volume, + annotation_stack=annotation_volume, + structures_list=structures_with_mesh, + meshes_dict=meshes_dict, + working_dir=bg_root_dir, + atlas_packager=PACKAGER, + hemispheres_stack=None, # it is symmetric + cleanup_files=False, + compress=True, + scale_meshes=True, + # resolution_mapping=[2, 1, 0], + ) + print("Done. Atlas generated at: ", output_filename) \ No newline at end of file From d981c91e5480c281bc458bada2ca320d21ce5f4f Mon Sep 17 00:00:00 2001 From: Carlo Castoldi Date: Mon, 25 Nov 2024 23:08:53 +0100 Subject: [PATCH 2/5] fix orientation and resolution --- .../atlas_scripts/kim_devccf_mouse.py | 15 ++++++--------- 1 file changed, 6 insertions(+), 9 deletions(-) diff --git a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py index baac91db..d4586645 100644 --- a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py +++ b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py @@ -22,7 +22,7 @@ AGE = "P14" ATLAS_LINK = "https://kimlab.io/brain-map/DevCCF/" CITATION = "Kronman, F.N., Liwang, J.K., Betty, R. et al. 2024, https://doi.org/10.1038/s41467-024-53254-w" -ORIENTATION = "posterior superior left" +ORIENTATION = ["left", "superior", "posterior"] ROOT_ID = 15564 RESOLUTION_UM = 20 VERSION = 1 @@ -130,7 +130,7 @@ def create_meshes(download_dir_path: str|Path, smooth = False # smooth meshes after creation start = time.time() if PARALLEL: - pool = mp.Pool(mp.cpu_count() - 2) + pool = mp.Pool(min(mp.cpu_count() - 2, 16)) try: pool.map( create_region_mesh, @@ -180,7 +180,6 @@ def create_meshes(download_dir_path: str|Path, ) return meshes_dir_path - def create_mesh_dict(structures, meshes_dir_path): meshes_dict = dict() structures_with_mesh = [] @@ -195,17 +194,14 @@ def create_mesh_dict(structures, meshes_dir_path): if mesh_path.stat().st_size < 512: print(f"obj file for {s} is too small, ignoring it.") continue - structures_with_mesh.append(s) meshes_dict[s["id"]] = mesh_path - print( f"In the end, {len(structures_with_mesh)} " "structures with mesh are kept" ) return meshes_dict, structures_with_mesh - if __name__ == "__main__": atlas_root = "/home/castoldi/Downloads/DevCCFv1" # Generated atlas path: @@ -213,11 +209,12 @@ def create_mesh_dict(structures, meshes_dir_path): bg_root_dir = DEFAULT_WORKDIR/ATLAS_NAME download_dir_path = bg_root_dir/"downloads" download_dir_path.mkdir(exist_ok=True, parents=True) - + structures = get_ontology(atlas_root) # save regions list json: with open(download_dir_path/"structures.json", "w") as f: json.dump(structures, f) + annotation_volume = get_annotations(atlas_root, AGE, RESOLUTION_UM) reference_volume = get_reference(atlas_root, AGE, RESOLUTION_UM, "LSFM") # Create meshes: @@ -228,7 +225,7 @@ def create_mesh_dict(structures, meshes_dir_path): annotation_volume, ROOT_ID ) - + # meshes_dir_path = download_dir_path/"meshes" meshes_dict, structures_with_mesh = create_mesh_dict( structures, meshes_dir_path ) @@ -240,7 +237,7 @@ def create_mesh_dict(structures, meshes_dir_path): citation=CITATION, atlas_link=ATLAS_LINK, species=SPECIES, - resolution=RESOLUTION_UM, + resolution=(RESOLUTION_UM,)*3, orientation=ORIENTATION, root_id=ROOT_ID, reference_stack=reference_volume, From 8b3783e168486cb834f31b528e9d7369b3b9edaa Mon Sep 17 00:00:00 2001 From: Carlo Castoldi Date: Fri, 20 Dec 2024 20:40:14 +0100 Subject: [PATCH 3/5] DevCCFv1 atlases are now fetched online and generated for all timestamps --- .../atlas_scripts/kim_devccf_mouse.py | 226 ++++++++---------- 1 file changed, 103 insertions(+), 123 deletions(-) diff --git a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py index d4586645..f07b7d12 100644 --- a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py +++ b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py @@ -2,24 +2,26 @@ import multiprocessing as mp import time from pathlib import Path +from os import listdir, path +import pooch import numpy as np import pandas as pd from brainglobe_utils.IO.image import load_nii from rich.progress import track +from brainglobe_atlasapi import utils from brainglobe_atlasapi.atlas_generation.mesh_utils import ( Region, create_region_mesh, ) from brainglobe_atlasapi.atlas_generation.wrapup import wrapup_atlas_from_data -# from brainglobe_atlasapi.config import DEFAULT_WORKDIR +from brainglobe_atlasapi.config import DEFAULT_WORKDIR from brainglobe_atlasapi.structure_tree_util import get_structures_tree ATLAS_NAME = "kim_devccf_mouse" # multiple versions of the same atlas SPECIES = "Mus musculus" -AGE = "P14" ATLAS_LINK = "https://kimlab.io/brain-map/DevCCF/" CITATION = "Kronman, F.N., Liwang, J.K., Betty, R. et al. 2024, https://doi.org/10.1038/s41467-024-53254-w" ORIENTATION = ["left", "superior", "posterior"] @@ -27,10 +29,7 @@ RESOLUTION_UM = 20 VERSION = 1 PACKAGER = "Carlo Castoldi " -ATLAS_FILE_URL = "https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/download.aspx?UniqueId=fe3d1692%2D94e4%2D4238%2Db6bc%2D95d18bcac022" -# curl 'https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/download.aspx?UniqueId=fe3d1692%2D94e4%2D4238%2Db6bc%2D95d18bcac022' -H 'User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:132.0) Gecko/20100101 Firefox/132.0' -H 'Accept: text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8' -H 'Accept-Language: en-GB,en;q=0.5' -H 'Accept-Encoding: gzip, deflate, br, zstd' -H 'Referer: https://pennstateoffice365-my.sharepoint.com/personal/yuk17_psu_edu/_layouts/15/onedrive.aspx?id=%2Fpersonal%2Fyuk17%5Fpsu%5Fedu%2FDocuments%2FProjects%2FDevelopmental%5FAtlasing%2FDevCCF%5FTeamFileshare%2FDevCCF%5FTemplates%2FDevCCFv1%2B%2FDevCCFv1' -H 'Upgrade-Insecure-Requests: 1' -H 'Sec-Fetch-Dest: iframe' -H 'Sec-Fetch-Mode: navigate' -H 'Sec-Fetch-Site: same-origin' -H 'Connection: keep-alive' -H 'Cookie: FedAuth=77u/PD94bWwgdmVyc2lvbj0iMS4wIiBlbmNvZGluZz0idXRmLTgiPz48U1A+VjEzLDBoLmZ8bWVtYmVyc2hpcHx1cm4lM2FzcG8lM2Fhbm9uIzc1Yjc1NjA0NGUxZDI5ZDU5NjFmMjZjODc0MzFlZmQwNmU4ZWM2NzYwMTJmZGEyNjgxYzg4Zjc5ODdmMGY2YTAsMCMuZnxtZW1iZXJzaGlwfHVybiUzYXNwbyUzYWFub24jNzViNzU2MDQ0ZTFkMjlkNTk2MWYyNmM4NzQzMWVmZDA2ZThlYzY3NjAxMmZkYTI2ODFjODhmNzk4N2YwZjZhMCwxMzM3NjQ5NTA5NzAwMDAwMDAsMCwxMzM3NjY3MzA4NzYxOTQyMTMsMC4wLjAuMCwyNTgsN2NmNDhkNDUtM2RkYi00Mzg5LWE5YzEtYzExNTUyNmViNTJlLCwsZWU2YTY1YTEtMDBmNS02MDAwLWZiMDYtYTgzYjk1Yjg3NTIzLDkxYzI2NWExLTMwNDctNjAwMC1mYjA2LWFiY2IwMTQ2MTZmNSxDR1EvcVlGcGJFeVlJTHAvM2NHb3pnLDAsMCwwLCwsLDI2NTA0Njc3NDM5OTk5OTk5OTksMCwsLCwsLCwwLCwxOTI2MDksdVhlaFFKUGxlVmpOQ2Jha1VoR0Q2SXlGUVFrLFkzdXYyY1dWdUZMWElneHRSOW4vMVVValJ1UUNMd2wybzdZbE51S0RtNTBxOE16ODJVaCt6K0ZTL0NqUU0rNlZrc090RlkwY211OUJNWHllc2U3Z3dST1VMcXRDWEorV3M1NkVzcEZrdk1qTEI0VzNhbDVZeVRlZ0greURRaUFRbG8rVG82aUVmVndVODU5Z2RWME5GSUprenJSVUlCeWVDUThiMmhMRjhQNmlJTDZHYldBNXVwRnJpcTl6cml6THhheFBRUWtVRzRoUmJjMEgzaVMzMnVKekw1ckFuWnVZSFk1N2J4K2xZOGE5MUNVQ3BSQXJ3aU45NHZ5aXdaQWdlcDRXUkxwWlVWOU5NNFZwZ0xUN05zVUh2K0JBcGwxbnViY0FDVENlWElhbWNwclZ0eTBFdElpajVSeXR3UGk0OFRlSVlIQXVrRjhDUFRFVEdDcTlGQT09PC9TUD4=; FeatureOverrides_experiments=[]; MicrosoftApplicationsTelemetryDeviceId=8293cff3-d815-447c-b343-546e678a0450; MSFPC=GUID=79c2671112584c09b960ec538dca634d&HASH=79c2&LV=202411&V=4&LU=1732022021463; ai_session=O/d7EbWPiaNaTZQh0qUYpP|1732113248412|1732113248412' --output P14.zip -PARALLEL = False # disable parallel mesh extraction for easier debugging - +ATLAS_FILE_URL = "https://doi.org/10.6084/m9.figshare.26377171.v1" TIMEPOINTS = ( "E11.5", @@ -51,35 +50,38 @@ "MRI-T2" # MRI T2-weighted ) -def get_annotations(devccfv1_root: str|Path, - age: str, - resolution: int): - if not isinstance(devccfv1_root, Path): - devccfv1_root = Path(devccfv1_root) +def pooch_init(download_dir_path: Path): + utils.check_internet_connection() + + hash = None + registry = {a+".zip": hash for a in TIMEPOINTS} + registry["DevCCFv1_OntologyStructure.xlsx"] = hash + return pooch.create( + path=download_dir_path, #/ATLAS_NAME, + base_url="doi:10.6084/m9.figshare.26377171.v1/", + registry=registry + ) + +def fetch_animal(pooch_: pooch.Pooch, age: str): assert age in TIMEPOINTS, f"Unknown age timepoint: '{age}'" - assert resolution in RESOLUTIONS, f"Unknown resolution in µm: '{resolution}'" - annotations_path = devccfv1_root/age/f"{age}_DevCCF_Annotations_{resolution}um.nii.gz" + archive = age+".zip" + members = [ + f"{age}/{age.replace('.','-')}_DevCCF_Annotations_20um.nii.gz", + f"{age}/{age.replace('.','-')}_LSFM_20um.nii.gz" + ] + annotations_path, reference_path = pooch_.fetch(archive, + progressbar=True, + processor=pooch.Unzip(extract_dir=".", members=members) + ) + # archive_path: Path = (pooch_.path/archive) + # archive_path.unlink() annotations = load_nii(annotations_path, as_array=True) - return annotations - -def get_reference(devccfv1_root: str|Path, - age: str, - resolution: int, - modality: str): - if not isinstance(devccfv1_root, Path): - devccfv1_root = Path(devccfv1_root) - assert age in TIMEPOINTS, f"Unknown age timepoint: '{age}'" - assert resolution in RESOLUTIONS, f"Unknown resolution: '{resolution}'" - assert modality in MODALITIES, f"Unknown modality: '{modality}'" - reference_path = devccfv1_root/age/f"{age}_{modality}_{resolution}um.nii.gz" reference = load_nii(reference_path, as_array=True) - return reference + return annotations, reference -def get_ontology(devccfv1_root: str|Path): - if not isinstance(devccfv1_root, Path): - devccfv1_root = Path(devccfv1_root) - DevCCFv1_path = devccfv1_root/"DevCCFv1_OntologyStructure.xlsx" - xl = pd.ExcelFile(DevCCFv1_path) +def fetch_ontology(pooch_: pooch.Pooch): + devccfv1_path = pooch_.fetch("DevCCFv1_OntologyStructure.xlsx", progressbar=True) + xl = pd.ExcelFile(devccfv1_path) # xl.sheet_names # it has two excel sheets # 'DevCCFv1_Ontology', 'README' df = xl.parse("DevCCFv1_Ontology", header=1) @@ -105,12 +107,11 @@ def get_ontology(devccfv1_root: str|Path): del structure["b"] return structures -def create_meshes(download_dir_path: str|Path, +def create_meshes(output_path: str|Path, structures, annotation_volume, root_id): - if not isinstance(download_dir_path, Path): - download_dir_path = Path(download_dir_path) - meshes_dir_path = download_dir_path / "meshes" - meshes_dir_path.mkdir(exist_ok=True) + if not isinstance(output_path, Path): + output_path = Path(output_path) + output_path.mkdir(exist_ok=True) tree = get_structures_tree(structures) @@ -129,56 +130,35 @@ def create_meshes(download_dir_path: str|Path, # What fraction of the original number of vertices is to be kept. smooth = False # smooth meshes after creation start = time.time() - if PARALLEL: - pool = mp.Pool(min(mp.cpu_count() - 2, 16)) - try: - pool.map( - create_region_mesh, - [ - ( - meshes_dir_path, - node, - tree, - labels, - annotation_volume, - root_id, - closing_n_iters, - decimate_fraction, - smooth, - ) - for node in tree.nodes.values() - ], - ) - except mp.pool.MaybeEncodingError: - pass - else: - for node in track( - tree.nodes.values(), - total=tree.size(), - description="Creating meshes", - ): - # root_node = tree.nodes[root_id] - create_region_mesh( - ( - meshes_dir_path, - # root_node, - node, - tree, - labels, - annotation_volume, - root_id, - closing_n_iters, - decimate_fraction, - smooth, - ) + for node in track( + tree.nodes.values(), + total=tree.size(), + # list(tree.nodes.values())[:5], + # total=5, + description="Creating meshes", + ): + # root_node = tree.nodes[root_id] + create_region_mesh( + ( + output_path, + # root_node, + node, + tree, + labels, + annotation_volume, + root_id, + closing_n_iters, + decimate_fraction, + smooth, ) + ) print( "Finished mesh extraction in: ", round((time.time() - start) / 60, 2), " minutes", ) - return meshes_dir_path + return output_path def create_mesh_dict(structures, meshes_dir_path): meshes_dict = dict() @@ -203,53 +183,53 @@ def create_mesh_dict(structures, meshes_dir_path): return meshes_dict, structures_with_mesh if __name__ == "__main__": - atlas_root = "/home/castoldi/Downloads/DevCCFv1" - # Generated atlas path: - DEFAULT_WORKDIR = Path.home()/"brainglobe_workingdir" bg_root_dir = DEFAULT_WORKDIR/ATLAS_NAME download_dir_path = bg_root_dir/"downloads" download_dir_path.mkdir(exist_ok=True, parents=True) - - structures = get_ontology(atlas_root) + pooch_ = pooch_init(download_dir_path) + structures = fetch_ontology(pooch_) # save regions list json: - with open(download_dir_path/"structures.json", "w") as f: + with open(bg_root_dir/"structures.json", "w") as f: json.dump(structures, f) - annotation_volume = get_annotations(atlas_root, AGE, RESOLUTION_UM) - reference_volume = get_reference(atlas_root, AGE, RESOLUTION_UM, "LSFM") - # Create meshes: - print(f"Saving atlas data at {download_dir_path}") - meshes_dir_path = create_meshes( - download_dir_path, - structures, - annotation_volume, - ROOT_ID - ) - # meshes_dir_path = download_dir_path/"meshes" - meshes_dict, structures_with_mesh = create_mesh_dict( - structures, meshes_dir_path - ) - # Wrap up, compress, and remove file: - print("Finalising atlas") - output_filename = wrapup_atlas_from_data( - atlas_name=ATLAS_NAME, - atlas_minor_version=VERSION, - citation=CITATION, - atlas_link=ATLAS_LINK, - species=SPECIES, - resolution=(RESOLUTION_UM,)*3, - orientation=ORIENTATION, - root_id=ROOT_ID, - reference_stack=reference_volume, - annotation_stack=annotation_volume, - structures_list=structures_with_mesh, - meshes_dict=meshes_dict, - working_dir=bg_root_dir, - atlas_packager=PACKAGER, - hemispheres_stack=None, # it is symmetric - cleanup_files=False, - compress=True, - scale_meshes=True, - # resolution_mapping=[2, 1, 0], - ) - print("Done. Atlas generated at: ", output_filename) \ No newline at end of file + for age in TIMEPOINTS: + atlas_name = f"{ATLAS_NAME}_{age.replace('.', '-')}" + annotation_volume, reference_volume = fetch_animal(pooch_, age) + atlas_dir = bg_root_dir/atlas_name + atlas_dir.mkdir(exist_ok=True) + print(f"Saving atlas data at {atlas_dir}") + # Create meshes: + meshes_dir_path = atlas_dir/"meshes" + create_meshes( + meshes_dir_path, + structures, + annotation_volume, + ROOT_ID + ) + meshes_dict, structures_with_mesh = create_mesh_dict( + structures, meshes_dir_path + ) + # Wrap up, compress, and remove file: + print("Finalising atlas") + output_filename = wrapup_atlas_from_data( + atlas_name=atlas_name, + atlas_minor_version=VERSION, + citation=CITATION, + atlas_link=ATLAS_LINK, + species=SPECIES, + resolution=(RESOLUTION_UM,)*3, + orientation=ORIENTATION, + root_id=ROOT_ID, + reference_stack=reference_volume, + annotation_stack=annotation_volume, + structures_list=structures_with_mesh, + meshes_dict=meshes_dict, + working_dir=atlas_dir, + atlas_packager=PACKAGER, + hemispheres_stack=None, # it is symmetric + cleanup_files=False, + compress=True, + scale_meshes=True, + # resolution_mapping=[2, 1, 0], + ) + print("Done. Atlas generated at: ", output_filename) \ No newline at end of file From d0addd940c2d85ad6539f88aa801c5367abc2256 Mon Sep 17 00:00:00 2001 From: Carlo Castoldi Date: Sat, 21 Dec 2024 11:30:11 +0100 Subject: [PATCH 4/5] skip mesh creation if done already --- .../atlas_generation/atlas_scripts/kim_devccf_mouse.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py index f07b7d12..4d575a5a 100644 --- a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py +++ b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py @@ -137,6 +137,10 @@ def create_meshes(output_path: str|Path, # total=5, description="Creating meshes", ): + output_file = output_path/f"{node.identifier}.obj" + if output_file.exists(): + # print(f"mesh already existing: {output_file.exists()} - {output_file}") + continue # root_node = tree.nodes[root_id] create_region_mesh( ( From 176eb0a962991c846ae10fd587342f82f6132d5e Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Sat, 21 Dec 2024 10:29:38 +0000 Subject: [PATCH 5/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../atlas_scripts/kim_devccf_mouse.py | 115 +++++++++--------- 1 file changed, 59 insertions(+), 56 deletions(-) diff --git a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py index 4d575a5a..2cddf67b 100644 --- a/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py +++ b/brainglobe_atlasapi/atlas_generation/atlas_scripts/kim_devccf_mouse.py @@ -1,12 +1,10 @@ import json -import multiprocessing as mp import time from pathlib import Path -from os import listdir, path -import pooch import numpy as np import pandas as pd +import pooch from brainglobe_utils.IO.image import load_nii from rich.progress import track @@ -16,7 +14,6 @@ create_region_mesh, ) from brainglobe_atlasapi.atlas_generation.wrapup import wrapup_atlas_from_data - from brainglobe_atlasapi.config import DEFAULT_WORKDIR from brainglobe_atlasapi.structure_tree_util import get_structures_tree @@ -31,84 +28,93 @@ PACKAGER = "Carlo Castoldi " ATLAS_FILE_URL = "https://doi.org/10.6084/m9.figshare.26377171.v1" -TIMEPOINTS = ( - "E11.5", - "E13.5", - "E15.5", - "E18.5", - "P04", - "P14", - "P56" -) +TIMEPOINTS = ("E11.5", "E13.5", "E15.5", "E18.5", "P04", "P14", "P56") RESOLUTIONS = (20, 50) MODALITIES = ( - "LSFM", # Light Sheet Fluorescence Microscopy + "LSFM", # Light Sheet Fluorescence Microscopy "MRI-adc", # MRI Apparent Diffusion Coefficient "MRI-dwi", # MRI Difusion Weighted Imaging - "MRI-fa", # MRI Fractional Anisotropy + "MRI-fa", # MRI Fractional Anisotropy "MRI-MTR", # MRI Magnetization Transfer Ratio - "MRI-T2" # MRI T2-weighted + "MRI-T2", # MRI T2-weighted ) + def pooch_init(download_dir_path: Path): utils.check_internet_connection() hash = None - registry = {a+".zip": hash for a in TIMEPOINTS} + registry = {a + ".zip": hash for a in TIMEPOINTS} registry["DevCCFv1_OntologyStructure.xlsx"] = hash return pooch.create( - path=download_dir_path, #/ATLAS_NAME, + path=download_dir_path, # /ATLAS_NAME, base_url="doi:10.6084/m9.figshare.26377171.v1/", - registry=registry + registry=registry, ) + def fetch_animal(pooch_: pooch.Pooch, age: str): assert age in TIMEPOINTS, f"Unknown age timepoint: '{age}'" - archive = age+".zip" + archive = age + ".zip" members = [ f"{age}/{age.replace('.','-')}_DevCCF_Annotations_20um.nii.gz", - f"{age}/{age.replace('.','-')}_LSFM_20um.nii.gz" + f"{age}/{age.replace('.','-')}_LSFM_20um.nii.gz", ] - annotations_path, reference_path = pooch_.fetch(archive, - progressbar=True, - processor=pooch.Unzip(extract_dir=".", members=members) - ) + annotations_path, reference_path = pooch_.fetch( + archive, + progressbar=True, + processor=pooch.Unzip(extract_dir=".", members=members), + ) # archive_path: Path = (pooch_.path/archive) # archive_path.unlink() annotations = load_nii(annotations_path, as_array=True) reference = load_nii(reference_path, as_array=True) return annotations, reference + def fetch_ontology(pooch_: pooch.Pooch): - devccfv1_path = pooch_.fetch("DevCCFv1_OntologyStructure.xlsx", progressbar=True) + devccfv1_path = pooch_.fetch( + "DevCCFv1_OntologyStructure.xlsx", progressbar=True + ) xl = pd.ExcelFile(devccfv1_path) # xl.sheet_names # it has two excel sheets - # 'DevCCFv1_Ontology', 'README' + # 'DevCCFv1_Ontology', 'README' df = xl.parse("DevCCFv1_Ontology", header=1) df = df[["Acronym", "ID16", "Name", "Structure ID Path16", "R", "G", "B"]] - df.rename(columns={ - "Acronym": "acronym", - "ID16": "id", - "Name": "name", - "Structure ID Path16": "structure_id_path", - "R": "r", - "G": "g", - "B": "b" - }, inplace=True) + df.rename( + columns={ + "Acronym": "acronym", + "ID16": "id", + "Name": "name", + "Structure ID Path16": "structure_id_path", + "R": "r", + "G": "g", + "B": "b", + }, + inplace=True, + ) structures = list(df.to_dict(orient="index").values()) for structure in structures: if structure["acronym"] == "mouse": structure["acronym"] = "root" structure_path = structure["structure_id_path"] - structure["structure_id_path"] = [int(id) for id in structure_path.strip("/").split("/")] - structure["rgb_triplet"] = [structure["r"], structure["g"], structure["b"]] + structure["structure_id_path"] = [ + int(id) for id in structure_path.strip("/").split("/") + ] + structure["rgb_triplet"] = [ + structure["r"], + structure["g"], + structure["b"], + ] del structure["r"] del structure["g"] del structure["b"] return structures -def create_meshes(output_path: str|Path, - structures, annotation_volume, root_id): + +def create_meshes( + output_path: str | Path, structures, annotation_volume, root_id +): if not isinstance(output_path, Path): output_path = Path(output_path) output_path.mkdir(exist_ok=True) @@ -126,7 +132,7 @@ def create_meshes(output_path: str|Path, # Mesh creation closing_n_iters = 2 - decimate_fraction = 0.2 # 0.04 + decimate_fraction = 0.2 # 0.04 # What fraction of the original number of vertices is to be kept. smooth = False # smooth meshes after creation start = time.time() @@ -137,7 +143,7 @@ def create_meshes(output_path: str|Path, # total=5, description="Creating meshes", ): - output_file = output_path/f"{node.identifier}.obj" + output_file = output_path / f"{node.identifier}.obj" if output_file.exists(): # print(f"mesh already existing: {output_file.exists()} - {output_file}") continue @@ -164,6 +170,7 @@ def create_meshes(output_path: str|Path, ) return output_path + def create_mesh_dict(structures, meshes_dir_path): meshes_dict = dict() structures_with_mesh = [] @@ -186,30 +193,26 @@ def create_mesh_dict(structures, meshes_dir_path): ) return meshes_dict, structures_with_mesh + if __name__ == "__main__": - bg_root_dir = DEFAULT_WORKDIR/ATLAS_NAME - download_dir_path = bg_root_dir/"downloads" + bg_root_dir = DEFAULT_WORKDIR / ATLAS_NAME + download_dir_path = bg_root_dir / "downloads" download_dir_path.mkdir(exist_ok=True, parents=True) pooch_ = pooch_init(download_dir_path) structures = fetch_ontology(pooch_) # save regions list json: - with open(bg_root_dir/"structures.json", "w") as f: + with open(bg_root_dir / "structures.json", "w") as f: json.dump(structures, f) for age in TIMEPOINTS: atlas_name = f"{ATLAS_NAME}_{age.replace('.', '-')}" annotation_volume, reference_volume = fetch_animal(pooch_, age) - atlas_dir = bg_root_dir/atlas_name + atlas_dir = bg_root_dir / atlas_name atlas_dir.mkdir(exist_ok=True) print(f"Saving atlas data at {atlas_dir}") # Create meshes: - meshes_dir_path = atlas_dir/"meshes" - create_meshes( - meshes_dir_path, - structures, - annotation_volume, - ROOT_ID - ) + meshes_dir_path = atlas_dir / "meshes" + create_meshes(meshes_dir_path, structures, annotation_volume, ROOT_ID) meshes_dict, structures_with_mesh = create_mesh_dict( structures, meshes_dir_path ) @@ -221,7 +224,7 @@ def create_mesh_dict(structures, meshes_dir_path): citation=CITATION, atlas_link=ATLAS_LINK, species=SPECIES, - resolution=(RESOLUTION_UM,)*3, + resolution=(RESOLUTION_UM,) * 3, orientation=ORIENTATION, root_id=ROOT_ID, reference_stack=reference_volume, @@ -230,10 +233,10 @@ def create_mesh_dict(structures, meshes_dir_path): meshes_dict=meshes_dict, working_dir=atlas_dir, atlas_packager=PACKAGER, - hemispheres_stack=None, # it is symmetric + hemispheres_stack=None, # it is symmetric cleanup_files=False, compress=True, scale_meshes=True, # resolution_mapping=[2, 1, 0], ) - print("Done. Atlas generated at: ", output_filename) \ No newline at end of file + print("Done. Atlas generated at: ", output_filename)