Find all integral coordinates given top left corner and bottom right corner of a rectangle
Since you commented that using numpy is OK for you, you can simplify the whole rectangle computation using its advanced slicing capabilities:

&gt;&gt;&gt; # Let's define a square array
... 
&gt;&gt;&gt; a = np.array(range(25)).reshape((5, 5))
&gt;&gt;&gt; a
array([[ 0,  1,  2,  3,  4],
       [ 5,  6,  7,  8,  9],
       [10, 11, 12, 13, 14],
       [15, 16, 17, 18, 19],
       [20, 21, 22, 23, 24]])
&gt;&gt;&gt; # Advanced slicing allow us to extract blocks at once
... 
&gt;&gt;&gt; a[2:5, 1:3]
array([[11, 12],
       [16, 17],
       [21, 22]])


You just need to encapsulate this behaviour in a class for ease of use:

import numpy as np


class BattleField:
    def __init__(self, size):
        self.field = np.zeros((size, size), dtype=int)

    def receive_missile(self, power, top, left, bottom, right):
        self.field[top-1:bottom, left-1:right] ^= power

    def __str__(self):
        return str(self.field)


if __name__ == '__main__':
    missiles = [
        [3,3,1,3,2],
        [2,2,1,2,2],
        [3,1,1,2,3],
    ]
    field = BattleField(3)
    for missile in missiles:
        field.receive_missile(*missile)
    print(field)




Other points to improve from your code includes:


using more descriptive names than A, B, or C;
using functions instead of code at top-level for ease of testing/reuse;
using the if __name__ == '__main__' guard to avoid running some code when importing the file;
using direct parameters values instead of packing them in a list and unpacking them in the function call (product(*[range(N), range(N)]) => product(range(N), range(N)));
using dict.get with its default value (coords[j] = coords.get(j, 0) ^ Pi) to simplify the code.