diff --git a/gprMax/exceptions.py b/gprMax/exceptions.py
index aaf4de17..75a34d9c 100644
--- a/gprMax/exceptions.py
+++ b/gprMax/exceptions.py
@@ -23,7 +23,7 @@ from colorama import Fore
 
 init()
 
-#sys.tracebacklimit = None
+sys.tracebacklimit = None
 
 
 class GeneralError(ValueError):
diff --git a/gprMax/grid.py b/gprMax/grid.py
index f8664624..78a8b8ab 100644
--- a/gprMax/grid.py
+++ b/gprMax/grid.py
@@ -89,14 +89,14 @@ class FDTDGrid(Grid):
 
         # Get information about host machine
         self.hostinfo = None
-        
+
         # CPU - OpenMP threads
         self.nthreads = 0
-        
+
         # GPU
         # Threads per block
         self.tpb = (256, 1, 1)
-        
+
         # GPU object
         self.gpu = None
 
@@ -238,8 +238,10 @@ def dispersion_analysis(G):
 
                 # Built-in waveform
                 else:
-                    time = np.linspace(0, 1, G.iterations)
-                    time *= (G.iterations * G.dt)
+                    # Time to analyse waveform - 3*pulse_width as using entire
+                    # time window can result in demanding FFT
+                    waveform.calculate_coefficients()
+                    time = np.arange(0, 3 * waveform.chi, G.dt)
                     waveformvalues = np.zeros(len(time))
                     timeiter = np.nditer(time, flags=['c_index'])
 
diff --git a/gprMax/waveforms.py b/gprMax/waveforms.py
index c89dfc46..a380ec40 100644
--- a/gprMax/waveforms.py
+++ b/gprMax/waveforms.py
@@ -25,7 +25,7 @@ class Waveform(object):
     """Definitions of waveform shapes that can be used with sources."""
 
     types = ['gaussian', 'gaussiandot', 'gaussiandotnorm', 'gaussiandotdot', 'gaussiandotdotnorm', 'gaussianprime', 'gaussiandoubleprime', 'ricker', 'sine', 'contsine', 'impulse', 'user']
-    
+
     # Information about specific waveforms:
     #
     # gaussianprime and gaussiandoubleprime waveforms are the first derivative and second derivative of the 'base' gaussian
@@ -40,6 +40,19 @@ class Waveform(object):
         self.amp = 1
         self.freq = None
         self.uservalues = None
+        self.chi = 0
+        self.zeta = 0
+        self.delay = 0
+
+    def calculate_coefficients(self):
+        """Calculates coefficients (used to calculate values) for specific waveforms."""
+
+        if self.type == 'gaussian' or self.type == 'gaussiandot' or self.type == 'gaussiandotnorm' or self.type == 'gaussianprime' or self.type == 'gaussiandoubleprime':
+            self.chi = 1 / self.freq
+            self.zeta = 2 * np.pi**2 * self.freq**2
+        elif self.type == 'gaussiandotdot' or self.type == 'gaussiandotdotnorm' or self.type == 'ricker':
+            self.chi = np.sqrt(2) / self.freq
+            self.zeta = np.pi**2 * self.freq**2
 
     def calculate_value(self, time, dt):
         """Calculates value of the waveform at a specific time.
@@ -52,37 +65,35 @@ class Waveform(object):
             ampvalue (float): Calculated value for waveform.
         """
 
-        # Coefficients for specific waveforms
-        if self.type == 'gaussian' or self.type == 'gaussiandot' or self.type == 'gaussiandotnorm' or self.type == 'gaussianprime' or self.type == 'gaussiandoubleprime':
-            chi = 1 / self.freq
-            zeta = 2 * np.pi**2 * self.freq**2
-            delay = time - chi
-        elif self.type == 'gaussiandotdot' or self.type == 'gaussiandotdotnorm' or self.type == 'ricker':
-            chi = np.sqrt(2) / self.freq
-            zeta = np.pi**2 * self.freq**2
-            delay = time - chi
-
+        self.calculate_coefficients()
+        
         # Waveforms
         if self.type == 'gaussian':
-            ampvalue = np.exp(-zeta * delay**2)
+            delay = time - self.chi
+            ampvalue = np.exp(-self.zeta * delay**2)
 
         elif self.type == 'gaussiandot' or self.type == 'gaussianprime':
-            ampvalue = -2 * zeta * delay * np.exp(-zeta * delay**2)
+            delay = time - self.chi
+            ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2)
 
         elif self.type == 'gaussiandotnorm':
-            normalise = np.sqrt(np.exp(1) / (2 * zeta))
-            ampvalue = -2 * zeta * delay * np.exp(-zeta * delay**2) * normalise
+            delay = time - self.chi
+            normalise = np.sqrt(np.exp(1) / (2 * self.zeta))
+            ampvalue = -2 * self.zeta * delay * np.exp(-self.zeta * delay**2) * normalise
 
         elif self.type == 'gaussiandotdot' or self.type == 'gaussiandoubleprime':
-            ampvalue = 2 * zeta * (2 * zeta * delay**2 - 1) * np.exp(-zeta * delay**2)
+            delay = time - self.chi
+            ampvalue = 2 * self.zeta * (2 * self.zeta * delay**2 - 1) * np.exp(-self.zeta * delay**2)
 
         elif self.type == 'gaussiandotdotnorm':
-            normalise = 1 / (2 * zeta)
-            ampvalue = 2 * zeta * (2 * zeta * delay**2 - 1) * np.exp(-zeta * delay**2) * normalise
+            delay = time - self.chi
+            normalise = 1 / (2 * self.zeta)
+            ampvalue = 2 * self.zeta * (2 * self.zeta * delay**2 - 1) * np.exp(-self.zeta * delay**2) * normalise
 
         elif self.type == 'ricker':
-            normalise = 1 / (2 * zeta)
-            ampvalue = - (2 * zeta * (2 * zeta * delay**2 - 1) * np.exp(-zeta * delay**2)) * normalise
+            delay = time - self.chi
+            normalise = 1 / (2 * self.zeta)
+            ampvalue = - (2 * self.zeta * (2 * self.zeta * delay**2 - 1) * np.exp(-self.zeta * delay**2)) * normalise
 
         elif self.type == 'sine':
             ampvalue = np.sin(2 * np.pi * self.freq * time)