compute the cell type specificity, defined as:
    ::
    
        T_yX = K_yX / K_YX
            = F_y * k_yX / sum_y(F_y*k_yX)
15120
def get_T_yX(fname, y, y_in_Y, x_in_X, F_y):
    '''
     
    
    '''
    def _get_k_yX_mul_F_y(y, y_index, X_index):
        # Load data from json dictionary
        f = open(fname, 'r')
        data = json.load(f)
        f.close()
    
        #init variables
        k_yX = 0.
        
        for l in [str(key) for key in data['data'][y]['syn_dict'].keys()]:
            for x in x_in_X[X_index]:
                p_yxL = data['data'][y]['syn_dict'][l][x] / 100.
                k_yL = data['data'][y]['syn_dict'][l]['number of synapses per neuron']
                k_yX +=  p_yxL * k_yL
        
        return k_yX * F_y[y_index]


    #container
    T_yX = np.zeros((len(y), len(x_in_X)))
    
    #iterate over postsynaptic cell types
    for i, y_value in enumerate(y):
        #iterate over presynapse population inds
        for j in range(len(x_in_X)):
            k_yX_mul_F_y = 0
            for k, yy in enumerate(sum(y_in_Y, [])):                
                if y_value in yy:
                    for yy_value in yy:
                        ii = np.where(np.array(y) == yy_value)[0][0]
                        k_yX_mul_F_y += _get_k_yX_mul_F_y(yy_value, ii, j)
            
            
            if k_yX_mul_F_y != 0:
                T_yX[i, j] = _get_k_yX_mul_F_y(y_value, i, j) / k_yX_mul_F_y
            
    return T_yX
896018
function subdivideLoop(cells, positions) {
  var edges       = nskel(cells, 1)
  var e_incidence = top.incidence(edges, cells)
  var dual        = top.dual(cells, positions.length)
  var npositions  = []
  var ncells      = []
  var e_indices   = new Array(edges.length)
  var v_indices   = new Array(positions.length)
  
  for(var i=0; i<e_indices.length; ++i) {
    e_indices[i] = -1
  }
  for(var i=0; i<v_indices.length; ++i) {
    v_indices[i] = -1
  }
  
  for(var f=0; f<cells.length; ++f) {
    var face = cells[f]
    
    for(var d=0; d<3; ++d) {
      var v = face[d]
      var u = face[(d+1)%3]
      e[0] = v
      e[1] = u
      var e_ptr = top.findCell(edges, e)
      var e_idx = e_indices[e_ptr]
      
      if(e_idx >= 0) {
        e_verts[d] = e_idx
      } else {
        //Compute edge-vertex
        var wing = e_incidence[e_ptr]
        var v0 = positions[u]
        var v1 = positions[v]
        var vertex = new Array(3)
        
        if(wing.length === 2) {
          var v2 = positions[opposite(u, v, cells[wing[0]])]
          var v3 = positions[opposite(u, v, cells[wing[1]])]
          
          for(var i=0; i<3; ++i) {
            vertex[i] = (3.0 * (v0[i] + v1[i]) + v2[i] + v3[i]) / 8.0
          }
        } else {
          for(var i=0; i<3; ++i) {
            vertex[i] = 0.5 * (v0[i] + v1[i])
          }
        }
        
        //Store vertex and continue
        e_indices[e_ptr] = e_verts[d] = npositions.length
        npositions.push(vertex)
      }
      
      if(v_indices[v] >= 0) {
        v_verts[d] = v_indices[v]
      } else {
        //Compute vertex-vertex weight
        
        //First, extract vertex neighborhood (slow and stupid here)
        var star = dual[v]
        var nbhd = [v]
        for(var i=0; i<star.length; ++i) {
          var tri = cells[star[i]]
          for(var j=0; j<3; ++j) {
            if(nbhd.indexOf(tri[j]) !== -1) {
              nbhd.push(tri[j])
            }
          }
        }
        
        //Next, compute weights
        var beta = (star.length === 3 ? 3.0/16.0 : 3.0/(8.0*star.length) )
        var center_weight = 1.0 - star.length * beta
        
        //Finally sum up weights
        var pos  = positions[v]
        var vertex = [0.1,0.1,0.1]
        for(var i=0; i<3; ++i) {
          vertex[i] = center_weight * pos[i]
        }
        for(var i=1; i<nbhd.length; ++i) {
          var p = positions[nbhd[i]]
          for(var j=0; j<3; ++j) {
            vertex[j] += beta * p[j]
          }
        }

        //Store result and continue        
        v_verts[d] = v_indices[v] = npositions.length
        npositions.push(vertex)
      }
    }
    
    //Add subdivided faces
    ncells.push([v_verts[0], e_verts[0], e_verts[2]],
                [e_verts[0], e_verts[1], e_verts[2]],
                [e_verts[0], v_verts[1], e_verts[1]],
                [e_verts[1], v_verts[2], e_verts[2]])
  }
  
  return { positions: npositions, cells: ncells }
}
868409
public static int resolveCellType(ICellValue c) {
        if (CellValue.BLANK == c) { return CELL_TYPE_BLANK; }

        if (c.getType() == Boolean.class) { return CELL_TYPE_BOOLEAN; }
        if (c.getType() == Double.class) { return CELL_TYPE_NUMERIC; }

        if (c.getType() == String.class) {
            String val = (String) c.get();

            if (val.startsWith(FORMULA_PREFIX)) { return CELL_TYPE_FORMULA; }
            if (isErrorValue(val)) { return CELL_TYPE_ERROR; }
            return CELL_TYPE_STRING;
        }

        throw new CalculationEngineException(String.format("CellVale's type %s is not supported.", c.get().getClass().getSimpleName()));
    }
186926
function dual(cells, vertex_count) {
  if(!vertex_count) {
    return incidence(unique(skeleton(cells, 0)), cells, 0)
  }
  var res = new Array(vertex_count)
  for(var i=0; i<vertex_count; ++i) {
    res[i] = []
  }
  for(var i=0, len=cells.length; i<len; ++i) {
    var c = cells[i]
    for(var j=0, cl=c.length; j<cl; ++j) {
      res[c[j]].push(i)
    }
  }
  return res
}
1235852
public CellType getCelltype() {
    if (CellTypeProteinConcentration_Type.featOkTst && ((CellTypeProteinConcentration_Type)jcasType).casFeat_celltype == null)
      jcasType.jcas.throwFeatMissing("celltype", "ch.epfl.bbp.uima.types.CellTypeProteinConcentration");
    return (CellType)(jcasType.ll_cas.ll_getFSForRef(jcasType.ll_cas.ll_getRefValue(addr, ((CellTypeProteinConcentration_Type)jcasType).casFeatCode_celltype)));}
439996
def cell(self, row_idx, col_idx):
        
        cell_idx = col_idx + (row_idx * self._column_count)
        return self._cells[cell_idx]
320331
def extractInputForTP(self, tm):
    

    # bursting cells in layer 4
    burstingColumns = tm.activeState["t"].sum(axis=1)
    burstingColumns[ burstingColumns < tm.cellsPerColumn ] = 0
    burstingColumns[ burstingColumns == tm.cellsPerColumn ] = 1
    # print "Bursting column indices=",burstingColumns.nonzero()[0]

    # correctly predicted cells in layer 4
    correctlyPredictedCells = numpy.zeros(self._inputDimensions).astype(realDType)
    idx = (tm.predictedState["t-1"] + tm.activeState["t"]) == 2
    idx = idx.reshape(self._inputDimensions)
    correctlyPredictedCells[idx] = 1.0
    # print "Predicted->active cell indices=",correctlyPredictedCells.nonzero()[0]

    # all currently active cells in layer 4
    spInputVector = tm.learnState["t"].reshape(self._inputDimensions)
    # spInputVector = tm.activeState["t"].reshape(self._inputDimensions)

    return (correctlyPredictedCells, spInputVector, burstingColumns)
17328
def _setCellClass (self):
        ''' '''
        from .. import sim
        
        # Check whether it's a NeuroML2 based cell
        if 'originalFormat' in self.tags:
            if self.tags['originalFormat'] == 'NeuroML2':
                self.cellModelClass = sim.NML2Cell
            if self.tags['originalFormat'] == 'NeuroML2_SpikeSource':
                self.cellModelClass = sim.NML2SpikeSource

        else:
            # set cell class: CompartCell for compartmental cells of PointCell for point neurons (NetStims, IntFire1,...)
            try: # check if cellModel corresponds to an existing point process mechanism; if so, use PointCell

                    tmp = getattr(h, self.tags['cellModel'])
                    self.cellModelClass = sim.PointCell
                    excludeTags = ['pop', 'cellModel', 'cellType', 'numCells', 'density', 'cellsList',
                                'xRange', 'yRange', 'zRange', 'xnormRange', 'ynormRange', 'znormRange', 'vref', 'spkTimes']
                    params = {k: v for k,v in self.tags.items() if k not in excludeTags}
                    self.tags['params'] = params
                    for k in self.tags['params']: self.tags.pop(k)
                    sim.net.params.popTagsCopiedToCells.append('params')
            except:
                if getattr(self.tags, 'cellModel', None) in ['NetStim', 'VecStim', 'IntFire1', 'IntFire2', 'IntFire4']:
                    print('Warning: could not find %s point process mechanism required for population %s' % (self.tags['cellModel'], self.tags['pop']))
                self.cellModelClass = sim.CompartCell
648951
def hkm_fc(fdata, Nmax, m, s):
    

    f = fdata[:, m]
    L1 = f.size
    MM = int(L1 / 2)
    Q = s.size

    ff = np.zeros(Q, dtype=np.complex128)
    for n in xrange(MM, L1):
        ff[n] = f[n - MM]

    for n in xrange(0, MM):
        ff[n] = f[n + MM]

    # For larger problems, this speeds things up pretty good.
    F = np.fft.fft(ff)
    S = np.fft.fft(s)
    out = 4 * np.pi * np.fft.ifft(F * S)

    return out[0:Nmax + 1]
1209194
def _cell_scalar(self, name=None):
        
        if name is None:
            # use active scalar array
            field, name = self.active_scalar_info
            if field != CELL_DATA_FIELD:
                raise RuntimeError('Must specify an array to fetch.')

        vtkarr = self.GetCellData().GetArray(name)
        if vtkarr is None:
            raise AssertionError('({}) is not a cell scalar'.format(name))

        # numpy does not support bit array data types
        if isinstance(vtkarr, vtk.vtkBitArray):
            vtkarr = vtk_bit_array_to_char(vtkarr)
            if name not in self._cell_bool_array_names:
                self._cell_bool_array_names.append(name)

        array = vtk_to_numpy(vtkarr)
        if array.dtype == np.uint8 and name in self._cell_bool_array_names:
            array = array.view(np.bool)
        return array
1059003
def _process_cell(i, state, finite=False):
    
    op_1 = state[i - 1]
    op_2 = state[i]
    if i == len(state) - 1:
        if finite:
            op_3 = state[0]
        else:
            op_3 = 0
    else:
        op_3 = state[i + 1]
    result = 0
    for i, val in enumerate([op_3, op_2, op_1]):
        if val:
            result += 2**i
    return result
1167739
def celltype(row, col, sheet = nil)
      sheet ||= default_sheet
      read_cells(sheet)
      row, col = normalize(row, col)
      @formula[sheet][[row, col]] ? :formula : @cell_type[sheet][[row, col]]
    end
469458
def dr( self, cell_lengths ):
        
        half_cell_lengths = cell_lengths / 2.0
        this_dr = self.final_site.r - self.initial_site.r
        for i in range( 3 ):
            if this_dr[ i ] > half_cell_lengths[ i ]:
                this_dr[ i ] -= cell_lengths[ i ]
            if this_dr[ i ] < -half_cell_lengths[ i ]:
                this_dr[ i ] += cell_lengths[ i ]
        return this_dr
829475
def cells(self):
        

        n = 0

        for (order, cells) in self:
            n += len(cells)

        return n
115116
def noise3(self, x, y, z):
		
		# Skew the input space to determine which simplex cell we're in
		s = (x + y + z) * _F3
		i = floor(x + s)
		j = floor(y + s)
		k = floor(z + s)
		t = (i + j + k) * _G3
		x0 = x - (i - t) # "Unskewed" distances from cell origin
		y0 = y - (j - t)
		z0 = z - (k - t)

		# For the 3D case, the simplex shape is a slightly irregular tetrahedron. 
		# Determine which simplex we are in. 
		if x0 >= y0:
			if y0 >= z0:
				i1 = 1; j1 = 0; k1 = 0
				i2 = 1; j2 = 1; k2 = 0
			elif x0 >= z0:
				i1 = 1; j1 = 0; k1 = 0
				i2 = 1; j2 = 0; k2 = 1
			else:
				i1 = 0; j1 = 0; k1 = 1
				i2 = 1; j2 = 0; k2 = 1
		else: # x0 < y0
			if y0 < z0:
				i1 = 0; j1 = 0; k1 = 1
				i2 = 0; j2 = 1; k2 = 1
			elif x0 < z0:
				i1 = 0; j1 = 1; k1 = 0
				i2 = 0; j2 = 1; k2 = 1
			else:
				i1 = 0; j1 = 1; k1 = 0
				i2 = 1; j2 = 1; k2 = 0
		
		# Offsets for remaining corners
		x1 = x0 - i1 + _G3
		y1 = y0 - j1 + _G3
		z1 = z0 - k1 + _G3
		x2 = x0 - i2 + 2.0 * _G3
		y2 = y0 - j2 + 2.0 * _G3
		z2 = z0 - k2 + 2.0 * _G3
		x3 = x0 - 1.0 + 3.0 * _G3
		y3 = y0 - 1.0 + 3.0 * _G3
		z3 = z0 - 1.0 + 3.0 * _G3

		# Calculate the hashed gradient indices of the four simplex corners
		perm = self.permutation
		ii = int(i) % self.period
		jj = int(j) % self.period
		kk = int(k) % self.period
		gi0 = perm[ii + perm[jj + perm[kk]]] % 12
		gi1 = perm[ii + i1 + perm[jj + j1 + perm[kk + k1]]] % 12
		gi2 = perm[ii + i2 + perm[jj + j2 + perm[kk + k2]]] % 12
		gi3 = perm[ii + 1 + perm[jj + 1 + perm[kk + 1]]] % 12

		# Calculate the contribution from the four corners
		noise = 0.0
		tt = 0.6 - x0**2 - y0**2 - z0**2
		if tt > 0:
			g = _GRAD3[gi0]
			noise = tt**4 * (g[0] * x0 + g[1] * y0 + g[2] * z0)
		else:
			noise = 0.0
		
		tt = 0.6 - x1**2 - y1**2 - z1**2
		if tt > 0:
			g = _GRAD3[gi1]
			noise += tt**4 * (g[0] * x1 + g[1] * y1 + g[2] * z1)
		
		tt = 0.6 - x2**2 - y2**2 - z2**2
		if tt > 0:
			g = _GRAD3[gi2]
			noise += tt**4 * (g[0] * x2 + g[1] * y2 + g[2] * z2)
		
		tt = 0.6 - x3**2 - y3**2 - z3**2
		if tt > 0:
			g = _GRAD3[gi3]
			noise += tt**4 * (g[0] * x3 + g[1] * y3 + g[2] * z3)
		
		return noise * 32.0
382930
def find_primitive(cell, symprec=1e-5):
    
    lattice, positions, numbers = spg.find_primitive(cell.totuple(), symprec)
    if lattice is None:
        return None
    else:
        return Atoms(numbers=numbers,
                     scaled_positions=positions,
                     cell=lattice,
                     pbc=True)
170950
def _K_computations(self, X, X2=None):
        """"""
        self._lengthscales=self.mapping.f(X)
        self._lengthscales2=np.square(self._lengthscales)
        if X2==None:
            self._lengthscales_two = self._lengthscales
            self._lengthscales_two2 = self._lengthscales2
            Xsquare = np.square(X).sum(1)
            self._K_dist2 = -2.*tdot(X) + Xsquare[:, None] + Xsquare[None, :]
        else:
            self._lengthscales_two = self.mapping.f(X2)
            self._lengthscales_two2 = np.square(self._lengthscales_two)
            self._K_dist2 = -2.*np.dot(X, X2.T) + np.square(X).sum(1)[:, None] + np.square(X2).sum(1)[None, :]
        self._w2 = self._lengthscales2 + self._lengthscales_two2.T
        prod_length = self._lengthscales*self._lengthscales_two.T
        self._K_exponential = np.exp(-self._K_dist2/self._w2)
        self._K_dvar = np.sign(prod_length)*(2*np.abs(prod_length)/self._w2)**(self.input_dim/2.)*np.exp(-self._K_dist2/self._w2)