diff --git a/cscripts/csphagen.py b/cscripts/csphagen.py
index ac36f80..aaecf00 100755
--- a/cscripts/csphagen.py
+++ b/cscripts/csphagen.py
@@ -3,7 +3,7 @@
 # Computes the PHA spectra for source/background and ARF/RMF files using the
 # reflected region method
 #
-# Copyright (C) 2017-2019 Luigi Tibaldo
+# Copyright (C) 2017-2020 Luigi Tibaldo
 #
 # This program is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -374,7 +374,7 @@ class csphagen(ctools.csobservation):
                 # Determine number of background regions to skip
                 N_skip  = self['bkgregskip'].integer()
                 N_start = 1 + N_skip
-                N_lim   = 1 + 2*N_skip
+                N_lim   = 1 + 2 * N_skip
 
                 # Compute the angular separation of reflected regions wrt
                 # camera center. The factor 1.05 ensures background regions
@@ -397,26 +397,36 @@ class csphagen(ctools.csobservation):
                 # regions
                 else:
 
-                    # Log appending of reflected regions
-                    msg = ' Use %d reflected regions.' % (N-N_lim)
-                    self._log_string(gammalib.NORMAL, msg)
-
                     # Append reflected regions
-                    alpha = 360.0 / N
-                    for s in range(N_start, N - N_skip):
-                        dphi    = s * alpha
+                    alpha    = 360.0 / N
+                    dphi_max = 360.0 - alpha * (1 + N_skip)
+                    dphi     = alpha         * (1 + N_skip)
+                    while dphi <= dphi_max:
                         ctr_dir = pnt_dir.clone()
                         ctr_dir.rotate_deg(posang + dphi, offset)
                         region = gammalib.GSkyRegionCircle(ctr_dir, self._rad)
                         if self._has_exclusion:
                             if self._excl_reg.overlaps(region):
+
+                                # Signal region overlap
                                 msg = ' Reflected region overlaps with '\
                                       'exclusion region.'
                                 self._log_string(gammalib.EXPLICIT, msg)
+
+                                # If region overlaps with exclusion region
+                                # try to increment by 10% of angular step
+                                dphi += 0.1 * alpha
+
                             else:
                                 regions.append(region)
+                                dphi += alpha
                         else:
                             regions.append(region)
+                            dphi += alpha
+
+                    # Log number of reflected regions
+                    msg = ' Use %d reflected regions.' % (regions.size())
+                    self._log_string(gammalib.NORMAL, msg)
 
         # Return reflected regions
         return regions
@@ -670,6 +680,13 @@ class csphagen(ctools.csobservation):
                 # Set On/Off observation ID
                 onoff.id(obs.id())
 
+            # Otherwise log observation skipped
+            else:
+                msg = ' Skip because the number %d of regions ' \
+                      'for background estimation is smaller than ' \
+                      '"bkgregmin"=%d.' % (bkg_reg.size(), self['bkgregmin'].integer())
+                self._log_string(gammalib.NORMAL, msg)
+
         # Construct dictionary with results
         result = {'onoff'     : onoff,
                   'bkg_reg'   : bkg_reg,