test commit.

Author: Dilawar Singh

biophysics/CaConcBase.cpp

@@ -33,28 +33,28 @@ SrcFinfo1< double >* CaConcBase::concOut()
 
 const Cinfo* CaConcBase::initCinfo()
 {
-    ///////////////////////////////////////////////////////
     // Shared message definitions
-    ///////////////////////////////////////////////////////
-    static DestFinfo process( "process",
-                              "Handles process call",
-                              new ProcOpFunc< CaConcBase >( &CaConcBase::process ) );
-    static DestFinfo reinit( "reinit",
-                             "Handles reinit call",
-                             new ProcOpFunc< CaConcBase >( &CaConcBase::reinit ) );
+    static DestFinfo process("process",
+            "Handles process call",
+            new ProcOpFunc< CaConcBase >( &CaConcBase::process )
+            );
+
+    static DestFinfo reinit("reinit",
+            "Handles reinit call",
+            new ProcOpFunc< CaConcBase >( &CaConcBase::reinit )
+            );
 
     static Finfo* processShared[] =
     {
         &process, &reinit
     };
 
-    static SharedFinfo proc( "proc",
-                             "Shared message to receive Process message from scheduler",
-                             processShared, sizeof( processShared ) / sizeof( Finfo* ) );
+    static SharedFinfo proc("proc",
+            "Shared message to receive Process message from scheduler",
+            processShared, sizeof( processShared ) / sizeof( Finfo* )
+            );
 
-    ///////////////////////////////////////////////////////
     // Field definitions
-    ///////////////////////////////////////////////////////
     static ElementValueFinfo< CaConcBase, double > Ca( "Ca",
             "Calcium concentration.",
             &CaConcBase::setCa,

python/moose/neuroml2/reader.py

@@ -22,6 +22,7 @@
 
 import logging
 logger_ = logging.getLogger('moose.nml2')
+logger_.setLevel(logging.INFO)
 
 try:
     import neuroml.loaders as loaders
@@ -49,6 +50,9 @@ def _unique(ls):
             res.append(l)
     return res
 
+def _whichGate(chan):
+    return chan.name[-1]
+
 def _isConcDep(ct):
     """_isConcDep
     Check if componet is dependant on concentration. Most HHGates are
@@ -59,19 +63,22 @@ def _isConcDep(ct):
 
     :return: True if Component is depenant on conc, False otherwise.
     """
+    if not hasattr(ct, 'extends'):
+        return False
     if 'ConcDep' in ct.extends:
         return True
     return False
 
-def _findCaConcVariableName():
-    """_findCaConcVariableName
+
+def _findCaConc():
+    """_findCaConc
     Find a suitable CaConc for computing HHGate tables.
     This is a hack, though it is likely to work in most cases. 
     """
     caConcs = moose.wildcardFind( '/library/##[TYPE=CaConc]' )
-    assert len(caConcs) >= 1, "No moose.CaConc found. Currently moose \
+    assert len(caConcs) == 1, "No moose.CaConc found. Currently moose \
             supports HHChannel which depends only on moose.CaConc ."
-    return caConcs[0].name
+    return caConcs[0]
 
 
 def sarea(comp):
@@ -471,17 +478,36 @@ def calculateRateFn(self, ratefn, vmin, vmax, tablen=3000, vShift='0mV'):
             logger_.info("Using %s to evaluate rate" % ct.name)
             rate = []
             for v in tab:
-                # Note: MOOSE HHGate are either voltage of concentration
-                # dependant. Here we figure out if nml description of gate is
-                # concentration dependant or not.
-                if _isConcDep(ct):
-                    # Concentration dependant. Concentration can't be negative.
-                    # Find a suitable CaConc from the /library. Currently only 
-                    # Ca dependant channels are allowed.
-                    caConcName = _findCaConcVariableName()
-                    req_vars  = {caConcName:'%g'%max(0,v),'vShift':vShift,'temperature':self._getTemperature()}
-                else:
-                    req_vars.update( self._variables )
+                req_vars = {'v':'%sV'%v,'vShift':vShift,'temperature':self._getTemperature()}
+                req_vars.update(self._variables)
+                vals = pynml.evaluate_component(ct, req_variables=req_vars)
+                '''print(vals)'''
+                if 'x' in vals:
+                    rate.append(vals['x'])
+                if 't' in vals:
+                    rate.append(vals['t'])
+                if 'r' in vals:
+                    rate.append(vals['r'])
+            return np.array(rate)
+
+    def calculateRateFnCaDep(self, ratefn, ca, camin, camax, tablen=3000, vShift='0mV'):
+        """Returns A / B table from ngate.
+
+        FIXME: Merge with calculateRateFn
+        """
+        tab = np.linspace(camin, camax, tablen)
+        if self._is_standard_nml_rate(ratefn):
+            midpoint, rate, scale = map(
+                SI, (ratefn.midpoint, ratefn.rate, ratefn.scale))
+            return self.rate_fn_map[ratefn.type](tab, rate, scale, midpoint)
+
+        for ct in self.doc.ComponentType:
+            if ratefn.type != ct.name:
+                continue
+            logger_.info("Using %s to evaluate rate (caConc dependant)" % ct.name)
+            rate = []
+            for v in tab:
+                req_vars = {ca.name:'%sV'%v,'vShift':vShift,'temperature':self._getTemperature()}
                 vals = pynml.evaluate_component(ct, req_variables=req_vars)
                 '''print(vals)'''
                 if 'x' in vals:
@@ -591,19 +617,23 @@ def _is_standard_nml_rate(self, rate):
     def createHHChannel(self, chan, vmin=-150e-3, vmax=100e-3, vdivs=5000):
         mchan = moose.HHChannel('%s/%s' % (self.lib.path, chan.id))
         mgates = map(moose.element, (mchan.gateX, mchan.gateY, mchan.gateZ))
-        assert (len(chan.gate_hh_rates) <= 3
-                )  # We handle only up to 3 gates in HHCHannel
+
+        # We handle only up to 3 gates in HHCHannel
+        assert (len(chan.gate_hh_rates) <= 3)  
 
         if self.verbose:
             print('== Creating channel: %s (%s) -> %s (%s)' %
                   (chan.id, chan.gate_hh_rates, mchan, mgates))
         all_gates = chan.gates + chan.gate_hh_rates
         for ngate, mgate in zip(all_gates, mgates):
-            if mgate.name.endswith('X'):
+            if ngate is None:
+                continue
+            logger_.info('whichGate(mgate) %s %s' %(mgate, _whichGate(mgate)))
+            if _whichGate(mgate) == 'X':
                 mchan.Xpower = ngate.instances
-            elif mgate.name.endswith('Y'):
+            elif _whichGate(mgate) == 'Y':
                 mchan.Ypower = ngate.instances
-            elif mgate.name.endswith('Z'):
+            elif _whichGate(mgate) == 'Z':
                 mchan.Zpower = ngate.instance
             mgate.min = vmin
             mgate.max = vmax
@@ -638,14 +668,28 @@ def createHHChannel(self, chan, vmin=-150e-3, vmax=100e-3, vdivs=5000):
                         'Unknown Q10 scaling type %s: %s' %
                         (ngate.q10_settings.type, ngate.q10_settings))
 
-            if self.verbose:
-                print(
-                    ' === Gate: %s; %s; %s; %s; %s; scale=%s' %
+            logger_.info(' === Gate: %s; %s; %s; %s; %s; scale=%s' %
                     (ngate.id, mgate.path, mchan.Xpower, fwd, rev, q10_scale))
 
+            print(fwd, rev)
+            quit()
+
             if (fwd is not None) and (rev is not None):
-                alpha = self.calculateRateFn(fwd, vmin, vmax, vdivs)
-                beta = self.calculateRateFn(rev, vmin, vmax, vdivs)
+                # Note: MOOSE HHGate are either voltage of concentration
+                # dependant. Here we figure out if nml description of gate is
+                # concentration dependant or not. 
+                if not _isConcDep(fwd):
+                    alpha = self.calculateRateFn(fwd, vmin, vmax, vdivs)
+                else:
+                    ca = _findCaConc()
+                    alpha = self.calculateRateFnCaDep(fwd, ca, vmin, vmax, vdivs)
+
+                if not _isConcDep(rev):
+                    beta = self.calculateRateFnCaDep(rev, vmin, vmax, vdivs)
+                else:
+                    ca = _findCaConc()
+                    beta = self.calculateRateFn(rev, ca, vmin, vmax, vdivs)
+
                 mgate.tableA = q10_scale * (alpha)
                 mgate.tableB = q10_scale * (alpha + beta)
 
@@ -659,9 +703,8 @@ def createHHChannel(self, chan, vmin=-150e-3, vmax=100e-3, vdivs=5000):
                 mgate.tableA = q10_scale * (inf / tau)
                 mgate.tableB = q10_scale * (1 / tau)
 
-            if hasattr(ngate,
-                       'steady_state') and (ngate.time_course is None) and (
-                           ngate.steady_state is not None):
+            if hasattr(ngate, 'steady_state') and (ngate.time_course is None) \
+                    and (ngate.steady_state is not None):
                 inf = ngate.steady_state
                 tau = 1 / (alpha + beta)
                 if (inf is not None):