Index: sasmodels/convert.py
===================================================================
--- sasmodels/convert.py	(revision 0a4628d166f9d3d4ed9811ec8882eea0e9675a8e)
+++ sasmodels/convert.py	(revision 34d6cabda48bc12af0423cded9c76c0e4bf2572b)
@@ -10,4 +10,5 @@
     'broad_peak',
     'two_lorentzian',
+    "two_power_law",
     'gel_fit',
     'gauss_lorentz_gel',
Index: sasmodels/models/line.py
===================================================================
--- sasmodels/models/line.py	(revision 0a4628d166f9d3d4ed9811ec8882eea0e9675a8e)
+++ sasmodels/models/line.py	(revision 34d6cabda48bc12af0423cded9c76c0e4bf2572b)
@@ -8,5 +8,10 @@
 .. math::
 
-    I(q) = A + B*q
+    I(q) = A + B \cdot q
+
+.. note::
+    For 2D plots intensity has different definition than other shape independent models
+.. math::
+    I(q) = I(qx) \cdot I(qy)
 
 .. figure:: None
Index: sasmodels/models/two_power_law.py
===================================================================
--- sasmodels/models/two_power_law.py	(revision 8115d829fb43c9da36d80783fde0626955690f0e)
+++ sasmodels/models/two_power_law.py	(revision 34d6cabda48bc12af0423cded9c76c0e4bf2572b)
@@ -11,11 +11,11 @@
 
     I(q) = \begin{cases}
-    A \pow(qval,-m1) & q <= qc \\
-    C \pow(qval,-m2) & q > qc
+    A q^{-m1} + \text{background} & q <= qc \\
+    C q^{-m2} + \text{background} & q > qc
     \end{cases}
 
 where $qc$ = the location of the crossover from one slope to the other,
 $A$ = the scaling coefficent that sets the overall intensity of the lower Q power law region,
-$m1$ = power law exponent at low Q
+$m1$ = power law exponent at low Q,
 $m2$ = power law exponent at high Q
 
@@ -24,7 +24,7 @@
 
 .. math::
-    C = frac{A}{\pow(qc/-m1)\pow(qc/-m2)}
+    C = \frac{A}{qc^{-m1} qc^{-m2}}
 
-...note
+.. note::
     Be sure to enter the power law exponents as positive values!
 
@@ -50,4 +50,5 @@
 from numpy import power
 from numpy import sqrt
+from numpy import inf
 
 name = "two_power_law"
@@ -68,14 +69,14 @@
 # pylint: disable=bad-whitespace, line-too-long
 #            ["name", "units", default, [lower, upper], "type", "description"],
-parameters = [["coef_A",  "",     1.0, [-inf, inf], "", "scaling coefficent in low Q region"],
-              ["qc", "1/Ang", 0.04, [0, inf], "", "crossover location"],
-              ["power_1",    "",    1.0, [0, inf], "", "power law exponent at low Q"],
-              ["power_2",  "",      4.0, [0, inf], "", "power law exponent at high Q"],
+parameters = [["coefficent_1",  "",         1.0, [-inf, inf], "", "coefficent A in low Q region"],
+              ["crossover",     "1/Ang",    0.04,[0, inf], "", "crossover location"],
+              ["power_1",       "",         1.0, [0, inf], "", "power law exponent at low Q"],
+              ["power_2",       "",         4.0, [0, inf], "", "power law exponent at high Q"],
               ]
 # pylint: enable=bad-whitespace, line-too-long
 
 def Iq(q,
-       coef_A=1.0,
-       qc=0.04,
+       coefficent_1=1.0,
+       crossover=0.04,
        power_1=1.0,
        power_2=4.0,
@@ -84,22 +85,21 @@
     """
     :param q:                   Input q-value (float or [float, float])
-    :param coef_A:              Scaling coefficent at low Q
-    :param qc:                  Crossover location
+    :param coefA:               Scaling coefficent at low Q
+    :param crossover:           Crossover location
     :param power_1:             Exponent of power law function at low Q
     :param power_2:             Exponent of power law function at high Q
     :return:                    Calculated intensity
     """
+# pylint: disable=bad-whitespace
 
-# pylint: disable=bad-whitespace
-    if(g<=qc):
-        intensity = coef_A*power(q,-1.0*power_1)
+    if q<=crossover:
+        intensity = coefficent_1*power(q,-1.0*power_1)
     else:
-        coef_C = coef_A*power(qc,-1.0*power_1)/power(qc,-1.0*power_2)
-        intensity = coef_C*power(q,-1.0*power_2)
+        coefficent_2 = coefficent_1*power(crossover,-1.0*power_1)/power(crossover,-1.0*power_2)
+        intensity = coefficent_2*power(q,-1.0*power_2)
 
     return intensity
 
-Iq.vectorized = True  # Iq accepts an array of q values
-
+Iq.vectorized = False  # Iq accepts an array of q values
 
 def Iqxy(qx, qy, *args):
@@ -113,41 +113,41 @@
     return Iq(sqrt(qx**2 + qy**2), *args)
 
-Iqxy.vectorized = True  # Iqxy accepts an array of qx, qy values
+Iqxy.vectorized = False  # Iqxy accepts an array of qx, qy values
 
 demo = dict(scale=1, background=0.1,
-            coef_A=1,
-            qc=0.04,
+            coefficent_1=1.0,
+            crossover=0.04,
             power_1=1.0,
             power_2=4.0)
 
 oldname = "TwoPowerLawModel"
-oldpars = dict(background='background',
-                coef_A='coef_A',
-                qc='qc',
+oldpars = dict(coefficent_1='coef_A',
+                crossover='qc',
                 power_1='power1',
-                power_2='power2')
+                power_2='power2',
+                background='background')
 
 tests = [
 
     # Accuracy tests based on content in test/utest_extra_models.py
-    [{'coeff_A':    1.0,
-      'qc':         0.04,
+    [{'coefficent_1':     1.0,
+      'crossover':  0.04,
       'power_1':    1.0,
       'power_2':    4.0,
-      'background': 0.1,
+      'background': 0.0,
     }, 0.001, 1000],
 
-    [{'coeff_A':    1.0,
-      'qc':         0.04,
+    [{'coefficent_1':     1.0,
+      'crossover':  0.04,
       'power_1':    1.0,
       'power_2':    4.0,
-      'background': 0.1,
+      'background': 0.0,
     }, 0.150141, 0.125945],
 
-    [{'coeff_A':    1.0,
-      'qc':         0.04,
+    [{'coeffcent_1':    1.0,
+      'crossover':  0.04,
       'power_1':    1.0,
       'power_2':    4.0,
-      'background': 0.1,
-    }, [0.442528,0.0], 0.00166884],
+      'background': 0.0,
+    }, 0.442528, 0.00166884],
 ]
