diff --git a/dist/libhreels-2.0.3.tar.gz b/dist/libhreels-2.0.3.tar.gz
deleted file mode 100644
index ab2e5be8eaf56532ae1112e811b5f152a08f171c..0000000000000000000000000000000000000000
Binary files a/dist/libhreels-2.0.3.tar.gz and /dev/null differ
diff --git a/dist/libhreels-2.0.3-py3-none-any.whl b/dist/libhreels-2.1.1-py3-none-any.whl
similarity index 54%
rename from dist/libhreels-2.0.3-py3-none-any.whl
rename to dist/libhreels-2.1.1-py3-none-any.whl
index 680aec3e0ad5df26cf143b1b5324eb7bbd1bf3bf..9cef8c4c28c6f223f222a7858ab2f3a792420be1 100644
Binary files a/dist/libhreels-2.0.3-py3-none-any.whl and b/dist/libhreels-2.1.1-py3-none-any.whl differ
diff --git a/dist/libhreels-2.1.1.tar.gz b/dist/libhreels-2.1.1.tar.gz
new file mode 100644
index 0000000000000000000000000000000000000000..6c0ce29960a60509097b8e4d356de15acc0979c2
Binary files /dev/null and b/dist/libhreels-2.1.1.tar.gz differ
diff --git a/libhreels/calcHREELS20.png b/libhreels/calcHREELS20.png
index d0332c96bfb10e992d8b4283e2bbb260c759c8d7..e893ca493faf568556151a0beae221a6a7e9d40b 100644
Binary files a/libhreels/calcHREELS20.png and b/libhreels/calcHREELS20.png differ
diff --git a/libhreels/calcHREELS20.py b/libhreels/calcHREELS20.py
index f01ab6ee01ff71ee2577fbb1120ac7a93523072d..78be6a5f3e70a74e84f0f3ae46aec6c4f11bb9ac 100755
--- a/libhreels/calcHREELS20.py
+++ b/libhreels/calcHREELS20.py
@@ -79,6 +79,8 @@ def addDrude(wLOPlasma, gLOPlasma, material, gTOPlasma='None'):
             newMaterial['gTO'] += [gTOPlasma]
             newMaterial['wLO'] += [wLOPlasma]
             newMaterial['gLO'] += [gLOPlasma]
+            if newMaterial.get('Q'):
+                newMaterial['Q'] += [0]
             return newMaterial
     except:
         print('Cannot add Drude to material',material)
@@ -102,7 +104,7 @@ class lambin:
         self.neps = self.layers
         # name_size = self.layers
         self.name = []; self.thick=[]; self.listNOsci=[]; self.epsinf =[]; Q = []
-        allTO=[]; allgTO=[];  allgLO=[]; nDrude=0
+        allTO=[]; allgTO=[];  allgLO=[]; nDrude=0; Qdummy = []
         name2 = []
         for layer in film:
             try:
@@ -122,14 +124,22 @@ class lambin:
             allgTO.extend(layer[0]['gTO'])
             allTO.extend(layer[0]['wLO'])
             allgTO.extend(layer[0]['gLO'])
-            Q.extend(2*len(layer[0]['wTO'])*[10.])
+            qList = layer[0].get('Q')
+            if qList:
+                Q.extend(layer[0]['Q'])
+                Q.extend(len(layer[0]['Q'])*[0.])
+            else:
+                Q.extend(2* len(layer[0]['wTO'])*[0.])
             self.listNOsci.append(nTO)
-        
+
         if len(allTO)!=sum(self.listNOsci) or len(allgTO)!=sum(self.listNOsci):
             print('Error in materials: ', layer[0])
+        if len(Q)!=sum(self.listNOsci) :
+            print('Error in materials (Check Q): ', layer[0])
         self.wOsc = np.array(allTO)
         self.gOsc = np.array(allgTO)
         self.osc = np.array([self.wOsc, np.array(Q), self.gOsc])
+        # print('[self.wOsc, np.array(Q), self.gOsc]: \n',self.osc)
         return
 
     def calcSurfaceLoss(self,x):
@@ -193,21 +203,56 @@ def myMain():
     import matplotlib.pyplot as plt
     import numpy as np
     import os
+    from libhreels import dielectrics20 as dielectrics
 
     x = np.linspace(-100.,1000,400)
-    material = {'eps': 4.,
-                'wTO': [0],  'gTO': [20], 'wLO': [598.7], 'gLO': [20]}
+    material = {'eps': 1.,
+                'wTO': [-200, -750],  
+                'gTO': [  12,   8], 
+                'wLO': [   1,   1],     # this parameter is irrelvant if wTO is negativ
+                'gLO': [   1,   1],     # this parameter is irrelvant if wTO is negativ
+                'Q'  : [  10,  15]}
+    material2 = {'eps': 1.,
+                'wTO': [-200, -750],  
+                'gTO': [  12,   8], 
+                'wLO': [   1,   1],     # this parameter is irrelvant if wTO is negativ
+                'gLO': [   1,   1],     # this parameter is irrelvant if wTO is negativ
+                'Q'  : [  0.1,  15]}
+    material3 = {'eps': 1.,
+                'wTO': [-200, -750,   0],  
+                'gTO': [  12,   8,   50], 
+                'wLO': [   1,   1, 2000],     # this parameter is irrelvant if wTO is negativ
+                'gLO': [   1,   1,   50],     # this parameter is irrelvant if wTO is negativ
+                'Q'  : [  10,  15,    0]}
 
     film1 = lambin(film=[[material,10000.]])
-    film1.temperature = 100
-    xs, spectrum = film1.calcHREELS(x,normalized=True,areanormalized=False)
+    film2 = lambin(film=[[material2,10000.]])
+    film3 = lambin(film=[[material3,10000.]])
+    eps3= film3.calcEps(x)[0]
+    eps2= film2.calcEps(x)[0]
+    eps = film1.calcEps(x)[0]
+    plt.plot(x,np.real(eps))
+    plt.plot(x,np.real(eps2))
+    plt.plot(x,np.real(eps3))
 
-    plt.plot(xs[:-1],spectrum[:-1], label='normalized=False')
-    plt.ylabel('Surface Loss')
+    ############# Comparison with dielectrics ################
+    osci1  = dielectrics.simpleOscillator(material3['wTO'][0], 
+                                          material3['Q'  ][0],
+                                      gTO=material3['gTO'][0])
+    osci2  = dielectrics.simpleOscillator(material3['wTO'][1], 
+                                          material3['Q'  ][1],
+                                      gTO=material3['gTO'][1])
+    drude = dielectrics.drude(material3['wLO'][2],material3['gTO'][2],material3['gLO'][2])
+    epsInfinity = material3['eps']
+    eps_dielectrics = epsInfinity * (osci1(x) + osci2(x)  + (1 + drude(x)))
+    plt.plot(x,np.real(eps_dielectrics),linestyle='dotted')
+
+
+
+    plt.ylabel(r'$Re(\epsilon)$')
     plt.xlabel('Energy Loss (cm$^{-1}$)')
-    plt.ylim(bottom=0)
-    plt.title('Bulk plasmon at 600 cm-1')
-    # plt.legend(title=r'Plasma frequency')
+    plt.xlim(left=5)
+    # plt.ylim(-6000,6000)
 
     plt.text(0.99, 0.01,os.path.basename(__file__), fontsize=10, ha='right', va='bottom', transform=plt.gcf().transFigure)
     output_filename = os.path.splitext(__file__)[0] + '.png'
@@ -215,6 +260,43 @@ def myMain():
 
     plt.show()
 
+    plt.plot(x,np.imag(eps))
+    plt.plot(x,np.imag(eps2))
+    plt.plot(x,np.imag(eps3))
+    plt.plot(x,np.imag(eps_dielectrics))
+    plt.xlim(left=0)
+    plt.ylim(-1500,1500)
+    plt.show()
+
+    # plt.plot(x,np.imag(dielectrics.sigma(eps,x)))
+    # plt.plot(x,np.imag(dielectrics.sigma(eps2,x)))
+    # plt.plot(x,np.imag(dielectrics.sigma(eps3,x)))
+    # plt.xlim(left=0)
+    # plt.ylim(-1500,1500)
+    # plt.show()
+
+    xs, spectrum = film3.calcHREELS(x,normalized=True,areanormalized=False)
+    plt.plot(xs[:-1],spectrum[:-1], label='normalized=False')
+    plt.show()
+
+    material4 = {'eps': 1.,
+                'wTO': [ 200, 750,   0],  
+                'gTO': [  12,   8,   50], 
+                'wLO': [ 600, 950, 2000],     
+                'gLO': [  10,  10,   50],     
+                'Q'  : [  10,  15,    0]}
+
+    # material4 = {'eps': 1.,
+    #             'wTO': [ 200, 750,   0],  
+    #             'gTO': [  12,   8,   50], 
+    #             'wLO': [ 300, 950, 2000], 
+    #             'gLO': [  10,  10,   50]}
+
+
+    film4 = lambin(film=[[material4,10000.]])
+    xs, spectrum = film4.calcHREELS(x,normalized=True,areanormalized=False)
+    plt.plot(xs[:-1],spectrum[:-1], label='normalized=False')
+    plt.show()
 
 if __name__ == '__main__':
 	myMain()
\ No newline at end of file
diff --git a/libhreels/dielectrics20.py b/libhreels/dielectrics20.py
index 584fd8cb929ce097d89ffa0b316861751c003b75..5729e57f3349c7f74dd0e8f077aa3b082caa2fa4 100755
--- a/libhreels/dielectrics20.py
+++ b/libhreels/dielectrics20.py
@@ -50,11 +50,13 @@ def sigma(eps,w):
     return (eps-1)*w/1j
 
 def int_sigma(eps,w):       
+    from scipy.integrate import cumulative_trapezoid
     '''Cumulatively integrates the optical conductivity for a given eps'''
     # Note that one value less is returned than the number of w values.
     # The following factor depends on the unit cell volume:
     fac = 0.67e-06
-    return fac*cumulative_trapezoid(np.real(sigma(eps,x)), x=w)
+    xx = w
+    return fac*cumulative_trapezoid(np.real(sigma(eps,w)), x=xx)
 
 def plotDielectrics(x,eps, title=" ", plot_show=True):
     '''This method plots a given (x,eps)-dataset as Re/Im(eps), SurfaceLoss and Reflectivity 
diff --git a/libhreels/materials20.json b/libhreels/materials20.json
index 2a656aca10d901855e92dc0ffeec47e466533036..c5658b26459b6c862cffd348bf53d9c36636705c 100755
--- a/libhreels/materials20.json
+++ b/libhreels/materials20.json
@@ -192,6 +192,15 @@
         "reference": "Irani'71, Ehrenreich'62, Furtak'75",
         "name": "Silver"
     },
+    "Ag-2": {
+        "eps": 5.0,
+        "wTO": [0],
+        "gTO": [161.0],
+        "wLO": [71800],
+        "gLO": [161.0],
+        "reference": "Yang'15",
+        "name": "Silver"
+    },
     "Pt": {
         "eps": 1.0,
         "wTO": [0],
@@ -236,6 +245,15 @@
         "gLO": [100],
         "reference": "...",
         "name": "vacuum"
+    },
+    "EuO": {
+        "eps": 4,
+        "wTO": [196 ],
+        "gTO": [ 30 ],
+        "wLO": [432 ],
+        "gLO": [ 30 ],
+        "reference": "Pradip'16",
+        "name": "EuO"
     }
 }
 
diff --git a/libhreels/mod_doeels.mod b/libhreels/mod_doeels.mod
new file mode 100644
index 0000000000000000000000000000000000000000..b3eaf47f6869cd2ff0f442a0d0914c8c68f63f7d
Binary files /dev/null and b/libhreels/mod_doeels.mod differ
diff --git a/libhreels/myEels20.cpython-38-x86_64-linux-gnu.so b/libhreels/myEels20.cpython-38-x86_64-linux-gnu.so
index dcb69b0d97b9e04ad9db972c13d186e86da59d08..05548bd84bfc3f7339bfe8cda10b25fbf3f55258 100755
Binary files a/libhreels/myEels20.cpython-38-x86_64-linux-gnu.so and b/libhreels/myEels20.cpython-38-x86_64-linux-gnu.so differ
diff --git a/libhreels/myEels20.f90 b/libhreels/myEels20.f90
index a91925b231af858fc363684c704c78055570e483..9a8bf821a1502b61ec8e1da5b2dd94b2353decc7 100755
--- a/libhreels/myEels20.f90
+++ b/libhreels/myEels20.f90
@@ -982,7 +982,7 @@ subroutine seteps(nos, osc_size, osc, epsinf, wn, nLayer, eps)
         deno =deno * (osc(1,j)**2 - wn * dcmplx( wn, osc(3,j) ) )
       
       else if (osc(1,j) < 0.) then! Negative TO mode means: _Additive_ Lorentz oscillator with Q
-        addeps = addeps + osc(1,j)**2 * osc(2,j) /dcmplx(osc(1,j)**2 + wn2, wn*osc(3,j))
+        addeps = addeps + osc(1,j)**2 * osc(2,j) /dcmplx(osc(1,j)**2 - wn2, -1*wn*osc(3,j))  ! Sign of imaginary part changed (WFW)
       
       else                      ! osc(1,j) = 0   -> it is a Drude term
         addeps = addeps - dcmplx(osc(1,j+m)**2, wn*(osc(3,j)-osc(3,j+m))) /dcmplx(wn2, wn*osc(3,j))
diff --git a/pyproject.toml b/pyproject.toml
index f9420b9fd81c66855e78eb96b3889319c575c073..eb505897dc619adf7acc9e3f85afbcf9f2c77a86 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "libhreels"
-version = "2.0.3"
+version = "2.1.1"
 description = "Handling, simulating, and plotting HREELS and Auger spectroscopy data"
 authors = ["Wolf Widdra <wolf.widdra@physik.uni-halle.de>"]
 include = ["*./libhreels/*"]