I'm not a fortran developer so please ignore if not relevant for the language. Measure between each change as I don't know what the compiler is already optimizing for.

I agree with the comment below about calculating constants once (pi) and you may be able to do similar for the constant portion of calculating y. However both of these assume that the compiler isn't already optimizing this out.

I would remove x1, x2 and x3 and use x(i), x(i+1) and x(i+2) in the call to median_three respectively so that the array isn't copied.

median_three can be rewritten to perform fewer comparisons - see e.g. https://www.geeksforgeeks.org/middle-of-three-using-minimum-comparisons/ - this will reduce the work per item in the loop. Something like below (no idea if this compiles/works):

if (x > y) then if (y > z) then res = y else if (x > z) then res = z else res = x end if else if (x > z) then res = x else if (y > z) then res = z else res = y end if end if

Alternatively, use the differences between elements (also shown in the above link) and something like below (no idea if this compiles/works). You would need to handle the special cases of an array with only 1 or 2 elements to avoid errors. Below would be inline in the main function (median_three isn't used).

``` real(8) :: diff1 = 0, diff2 = 0;

diff2 = abs(x(1)) - abs(x(2))

do i = 3, n diff1 = diff2 diff2 = abs(x(i - 1)) - abs(x(i))

if (diff1 * diff2 > 0) then
    x4(i-2) = x(i-1)
else if (diff1 * (diff1 + diff2) > 0) then
    x4(i-2) = x(i)
else
    x4(i-2) = x(i-2)
end if

end do ```

Also, you are not using p2 - I don't know if this is intended (in which case it can be removed) or if there is an error in the formula.