Hi Colin
Just at first glance... I don't think that this is too away from a deconvolution (or, a high pass convolution). You may try to generate a sinthetic image with a radial profile following that R^-2 decay, and use that image (normalize it so the total value is 1, Average*number_of_pixels).
As I said, this may also be done as a highpass filter. Since you are analysing a neighbourhood, you may calculare the gradient at a givent point, and create a function that rescales the image acording to that gradient. Other alternative may be to work directly on the gradient domain, and then use the Poisson solver Georg suggested to get the image. If you follow an iterative approach, you may emulate the 1/R^2 model with some success.