# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 2
#     language: python
#     name: python2
# ---

# +
import os
from allensdk.core.mouse_connectivity_cache import MouseConnectivityCache
import numpy as np
import h5py
import time
import nrrd
from voxnet.conn2d import map_to_surface

drive_path = os.path.join(os.getenv('HOME'), 'work/allen/data/sdk_new_100')

# When downloading 3D connectivity data volumes, what resolution do you want (in microns)?  
# Options are: 10, 25, 50, 100
resolution_um = 10

# Drop list criterion, in percent difference
volume_fraction = 20

# The manifest file is a simple JSON file that keeps track of all of
# the data that has already been downloaded onto the hard drives.
# If you supply a relative path, it is assumed to be relative to your
# current working directory.
manifest_file = os.path.join(drive_path, "manifest.json")

# Start processing data
mcc = MouseConnectivityCache(manifest_file = manifest_file,
                             resolution = resolution_um)
ontology = mcc.get_ontology()
# Injection structure of interest
isocortex = ontology['Isocortex']

# open up a pandas dataframe of all of the experiments
experiments = mcc.get_experiments(dataframe = True, 
                                  injection_structure_ids = [isocortex['id'].values[0]], 
                                  cre = False)
print "%d total experiments" % len(experiments)

view_paths_fn = os.path.join(os.getenv('HOME'), 'work/allen/data/TopView/top_view_paths_10.h5')
view_paths_file = h5py.File(view_paths_fn, 'r')
view_lut = view_paths_file['view lookup'][:]
view_paths = view_paths_file['paths'][:]
view_paths_file.close()

## Compute size of each path to convert path averages to sums
norm_lut = np.zeros(view_lut.shape, dtype=int)
ind = np.where(view_lut != -1)
ind = zip(ind[0], ind[1])
for curr_ind in ind:
    curr_path_id = view_lut[curr_ind]
    curr_path = view_paths[curr_path_id, :]
    norm_lut[curr_ind] = np.sum(curr_path > 0)

t0 = time.time()
expt_drop_list = []
full_vols = []
flat_vols = []
# -

eid = experiments.iloc[5].id
row = experiments.iloc[5]
print eid
print row

    data_dir = os.path.join(os.getenv('HOME'),
                            "work/allen/data/sdk_new_100/experiment_%d/" % eid)
    # get and remap injection data
    print "getting injection density"
    in_d, in_info = mcc.get_injection_density(eid)
    print "mapping to surface"
    in_d_s = map_to_surface(in_d, view_lut, view_paths, scale = resolution_um/10., fun=np.mean)
    flat_vol = np.nansum(in_d_s * norm_lut) * (10./1000.)**3
    flat_vols.append(flat_vol)
    full_vol = np.nansum(in_d) * (10./1000.)**3
    full_vols.append(full_vol)


expt_union = mcc.get_experiment_structure_unionizes(eid, hemisphere_ids = [2], is_injection = True,
                                                        structure_ids = [ontology['grey']['id'].values[0]])

flat_vols

full_vols

full_vol2 = float(expt_union['projection_volume'])
print expt_union
print full_vol2

print np.nansum(in_d_s * norm_lut)

in_d_sum = map_to_surface(in_d, view_lut, view_paths, scale = resolution_um/10., fun=np.sum)
np.nansum(in_d_sum) * (10./1000.)**3

print np.abs(flat_vol - full_vol)/full_vol*100
print np.abs(flat_vol - full_vol2)/full_vol2*100

# +
import matplotlib.pyplot as plt
# %matplotlib inline
tmp = norm_lut.astype(float)
tmp[tmp == 0] = np.nan

fig = plt.figure(figsize = (20,20))
ax = fig.add_subplot(121)
h = ax.imshow(in_d_s)
fig.colorbar(h)

a2 = fig.add_subplot(122)
h2 = a2.imshow(np.sum(in_d, axis=1))
fig.colorbar(h2)
# -

experiments



StarcoderdataJupyterScriptsDedupFiltered