R1 and R2 both include foreground and background regions where the foreground is in focus and the background is not. That is probably confusing your algorithms.

There is an interesting discussion on this subject on StackOverflow (see attached link), but I think for your image with its limited depth of field, you will need to do some segmentation before applying some of the suggested methods (or divide into (say) 10x10 blocks (or even overlapping blocks) and evaluate each one.

http://stackoverflow.com/questions/7765810/is-there-a-way-to-detect-if-an-image-is-blurry

With regions that are out of focus you tend to get more shallow gradients around the edges, as is the case with R1 and R2. Start by looking for the edges, but then also consider the approach to those edges. For example, if you look of at the Fourier transform of neighborhoods around the edges of R1 and R2 you will get a lot of mid-range frequencies. Around R3, on the other hand, you will get low frequencies and high-frequencies with some noise in the middle. By examining a number of these neighborhoods, you should be able to identify what types of frequencies to look for.

I think it will be very challenging to solve this problem using only one image. 

If you can obtain another image after tilting the camera for some angle, then you can estimate the depth using techniques from stereo imaging and use this in addition to the focus information that can be obtained from  the FFT.

http://uk.mathworks.com/help/vision/examples/depth-estimation-from-stereo-video.html?refresh=true

Well focused region means edges in that region are sharp. If the region is entirely smooth then in principle you cannot tell whether it's well focused or not without further knowledge of the image. 

However, if I were you, I would try the following on your image see what happens. Do a high pass filtering with various types of filters (different sizes/scales and difference orders of vanishing moments, etc.). Then each pixel now has a vector value with each entry being the result of filtering at that location using each filter. After this step, you may do an unsupervised classification on the pixels (2 classes): each pixel is now a data point in an M-dimensional space with M being the number of filters you use. See if you can clearly identify R3.

I will depart from what was said above. You could conceive a simple metric to yield how “edgy” a region is (as hinted by Cory and Bin) in combination with what Daniel suggested above (block processing).

The following MATLAB code would work for your supplied image but I am not sure of its generalization:

AA=edge(A,'canny'); % where A is your input image

fun = @(block_struct) nnz(block_struct.data) * ones(size(block_struct.data));

I = blockproc(AA,[64 64],fun);

figure; imshow(I,[]);

Dear Abbas,

What is blockproc?

I've got this message

Undefined function or method 'blockproc' for input arguments of type 'function_handle'.

Dear Mohamed,

It is a MATLAB built in function which you suppose to have if you have the Image Processing toolbox installed. I am using MATLAB R2014a. You can check the following link for a detailed explanation:

http://se.mathworks.com/help/images/ref/blockproc.html

N.B. In old versions you may find the same function being given a different name (e.g., blkproc in some older versions).

I guess it would be nice to hear from the question's poster (Hossein) on how did it go with his experiments or at least to enrich this thread with his feedback.

thanks you all researchers and experts. I applied most of points you mentioned above but, unfortunately none of them are general and always there are some exceptions. Also, it depends on sliding window size. some focusing measures with larger or smaller window size, find R3 as well focused region.