atan2() implemented in bc

Dave Taylor thought that it would be a good idea to calculate in Perl the distance between two longitude / latitude points on the Earth. He also thought it would be awesome if the same could be implemented in bc.

I replied that one could easily transform in Perl the cost function from the TSP page. Dave Taylor pointed out that unfortunately atan2() is not available in bc and atan2() is used in the cost function.

“So what?” I thought. All one needs is the definition of atan2(). It turns out that one needs to define in bc the sign function and abs() too.

define sgn(x) {
        if (x == 0) {
                return(0);
        } else if (x < 0) {
                return(-1);
        } else if (x > 0) {
                return(1);
        }
}

define abs(x) {
        if (x < 0) {
                return(-1 * x);
        } else {
                return(x);
        }
}

define atan2(y, x) {
        auto pi, fi;

        pi = 4 * a(1);

        if (y == 0) {
                if (x > 0) {
                        return(0);
                } else if (x == 0) {
                        print "undefined\n";
                        halt;
                } else if (x < 0) {
                        return(pi);
                }
        }

        fi = a(abs(y/x));

        if (x > 0) {
                return(fi * sgn(y));
        } else if (x == 0) {
                return(pi * sgn(y) / 2);
        } else if (x < 0) {
                return((pi - fi) * sgn(y));
        }
}

Put the function definitions in a file (say atan2.bc) and call it from the command line as:

$ bc -l atan2.bc

The above code was written in less than ten minutes of time. I am sure that one can come up with better implementations (not calculating pi (== 4 * atan(1)) and then dividing it by 2 is an optimization directly observable). All in all it was a fun exercise, plus while searching I got to read about CORDIC.