A computer-implemented method of radiometric correction of a group of partially overlapping images captured by an imaging device, the method comprising: using a processing circuitry for: for each given image in the group of images: determining a respective subset of overlapping images, each image having an area overlapping with the given image; determining a first set of sampled points in the given image, each sampled point having a respective value; determining a respective target value for each sampled point in the first set of sampled points, based on the respective value of the sampled point and values of homologous points in other images in the respective subset of overlapping images; determining a second set of sampled points in non-overlapping areas in the given image; determining a respective target value for each sampled point in the second set based on the respective value; executing an optimization procedure comprising: (i) for each sampled point i in the first set and the second set, calculating a respective difference DP i between the value of the sampled point and a respective target value, thereby generating a collection of difference values and a data space with the target values, wherein one dimension of the data space represents the sample point values and the other dimensions of the data space represent the respective coordinates of the sampled points; (ii) while a DP score based on the collection of DP values exceeds a threshold: fitting at least one optimized function (f) of the coordinates that approximates the sample point values in the data space in a manner that the DP score is reduced, and adding the at least one optimized function to a function-group, which is a collection of one or more optimized functions that were fitted for the given image; determining a correction function, based on the function-group; applying the correction function on the given image to thereby generate a respective corrected image. for each given image in the group of images: determining a respective subset of overlapping images, each image having an area overlapping with the given image; determining a first set of sampled points in the given image, each sampled point having a respective value; determining a respective target value for each sampled point in the first set of sampled points, based on the respective value of the sampled point and values of homologous points in other images in the respective subset of overlapping images; determining a second set of sampled points in non-overlapping areas in the given image; determining a respective target value for each sampled point in the second set based on the respective value; executing an optimization procedure comprising: (i) for each sampled point i in the first set and the second set, calculating a respective difference DP i between the value of the sampled point and a respective target value, thereby generating a collection of difference values and a data space with the target values, wherein one dimension of the data space represents the sample point values and the other dimensions of the data space represent the respective coordinates of the sampled points; (ii) while a DP score based on the collection of DP values exceeds a threshold: fitting at least one optimized function (f) of the coordinates that approximates the sample point values in the data space in a manner that the DP score is reduced, and adding the at least one optimized function to a function-group, which is a collection of one or more optimized functions that were fitted for the given image; (i) for each sampled point i in the first set and the second set, calculating a respective difference DP i between the value of the sampled point and a respective target value, thereby generating a collection of difference values and a data space with the target values, wherein one dimension of the data space represents the sample point values and the other dimensions of the data space represent the respective coordinates of the sampled points; (ii) while a DP score based on the collection of DP values exceeds a threshold: fitting at least one optimized function (f) of the coordinates that approximates the sample point values in the data space in a manner that the DP score is reduced, and adding the at least one optimized function to a function-group, which is a collection of one or more optimized functions that were fitted for the given image; determining a correction function, based on the function-group; applying the correction function on the given image to thereby generate a respective corrected image.