1

I have a raster/grid representing e.g. terrain aspect or flow direction, coded in geographic degrees (0-360°, clockwise from the north). I would like to convert it in arithmetic radians (0-2π, counterclockwise from the east).

I have no issue with converting degrees to radians in a linear fashion, but how can I handle the rotation and direction change? I am looking at achieving this programmatically (R/Python) or through usual GIS (Arc/Q), although I am favoring using QGIS' raster calculator.

4
  • 1
    Please decide, which software you're working with and add the appropriate tag, thank you.
    – Erik
    Mar 8, 2022 at 13:25
  • Edited per your request.
    – M2FKXY
    Mar 8, 2022 at 13:27
  • 1
    You can just multtipy the degree values by (pi/180). Just understand that the radians are clockwise from N, same as the degree angles. I don't see what you have to "handle" regarding rotation and direction.
    – Micha
    Mar 8, 2022 at 16:07
  • My needs specifically require my radians file to be counterclockwise from the east - simply converting degrees to radians is not sufficient.
    – M2FKXY
    Mar 9, 2022 at 15:32

2 Answers 2

4

The aspect tool in QGIS has an option to output your aspects in trigonometric angle. This would give you 0 as due east, 180 as west, 90 as south, and 270 as north. Then just multiply the result by pi/180 to convert to theta.

enter image description here

I see QGIS 3.22 now has an IF available in the raster calculator (but I do not have experience working with IF statements in the QGIS raster calculator). See this link for using If statements. If you do not have access to the original surface for input into a new Aspect calculation you could adapt this Python code to give you theta counter-clockwise from east.

def TrueNorthAzimuthToUnitCircleRadians(degrees):
    '''degrees in north 0 clockwise to unit circle radians'''
    import math
    if degrees > 0 and degrees < 180:
        degrees = degrees + 180
    elif degrees > 180 and degrees <=359:
        degrees = degrees - 180
    elif degrees ==  360:
        degrees = 0
    theta = degrees * math.pi/180.0
    return theta
0

Thanks GBG. I might have some something poorly with your code as I did not get the expected result, but looking at your answer to pretty much the same question but the other way around I was able to figure it out with minimal adaptation.

A one liner in QGIS' raster calculator would be:

2 * 3.14159265359 - (if("degrees" < 90, "degrees" + 270, "degrees" - 90) * 3.14159265359 / 180)

degrees being your raster@band.

Which roughly translated to pseudocode would be:

function degreesGeoCWtoRadiansArithmCCW(degreesGeoCW):
    degreesRotated = (degreesGeoCW - 90) % 360
    radians = degreesRotated * pi / 180
    radiansArithmCCW = 2 * pi - radians 

    return radiansArithmCCW

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.