Minecraft projectile aimbot

package it.tristana.arrows.helper;

import org.bukkit.util.Vector;

import static java.lang.Math.pow;
import static java.lang.Math.log;

public class ArrowMath {

    private ArrowMath() {}

    public static Vector toHit(Vector start, Vector target, Vector targetVelocity, double projectileSpeed, double drag, double gravity) {
        double px0 = start.getX();
        double py0 = start.getY();
        double pz0 = start.getZ();

        double pxt = target.getX();
        double pyt = target.getY();
        double pzt = target.getZ();

        double vx = targetVelocity.getX();
        double vy = targetVelocity.getY();
        double vz = targetVelocity.getZ();

        double J = pxt - px0;
        double B = pzt - pz0;
        double D = drag - 1;
        double X = gravity * drag;
        double F = D * gravity - vy * sq(D);
        double G = - (pyt - py0) * sq(D);
        double H = 2 * J * vx;
        double I = 2 * B * vz;
        double V = - 2 * sq(D) * drag;
        double Z = sq(D) * sq(drag);
        double Q = 2 * F * X;
        double M = 2 * X * G;
        double W = 2 * F * G;

        double C = sq(projectileSpeed) * Z * sq(drag) - sq(X) * sq(drag);
        double O = sq(projectileSpeed) * V * sq(drag) + 2 * sq(X) * sq(drag) + M * sq(drag) - 2 * sq(projectileSpeed) * Z * drag + 2 * sq(X) * drag;
        double U = sq(projectileSpeed) * sq(D) * sq(drag) - sq(X) * sq(drag) - sq(G) * sq(drag) - M * sq(drag) - 2 * sq(projectileSpeed) * V * drag - 4 * sq(X) * drag - 2 * M * drag + sq(projectileSpeed) * Z - sq(X) - sq(D) * sq(J) * Z - sq(D) * sq(B) * Z;
        double N = - 2 * sq(projectileSpeed) * sq(D) * drag + 2 * sq(X) * drag + 2 * sq(G) * drag + 2 * M * drag + sq(projectileSpeed) * V + 2 * sq(X) + M - sq(D) * sq(J) * V - sq(D) * sq(B) * V;
        double T = - sq(F) - pow(D, 4) * sq(vx) - pow(D, 4) * sq(vz);
        double E = - Q - W - pow(D, 4) * H - pow(D, 4) * I;
        double R = - sq(F) * sq(drag) - sq(D) * sq(vx) * Z - sq(D) * sq(vz) * Z;
        double P = - Q * sq(drag) - W * sq(drag) - 2 * Q * drag - sq(D) * H * Z - sq(D) * I * Z;
        double L = Q * sq(drag);
        double A = 2 * sq(F) * drag - sq(D) * sq(vx) * V - sq(D) * sq(vz) * V;
        double Y = 2 * W * drag + Q + 2 * Q * drag - sq(D) * H * V - sq(D) * I * V;
        double S = sq(projectileSpeed) * sq(D) - sq(X) - sq(G) - M - pow(D, 4) * (sq(J) + sq(B));

        double min = -1;
        double[] times = { 5, 10, 25, 50, 100, 250 };
        for (double test : times) {
            double time = getNewtonRaphsonApproximation(C, O, U, N, T, E, R, P, L, A, Y, S, drag, test, 1e-12, 1000);
            if (time > 0 && (time < min || min == -1)) {
                min = time;
            }
        }

        if (min < 0) {
            return null;
        }

        double x0 = getHorizontalComponent(J, D, vx, drag, min);
        double y0 = getVerticalComponent(D, X, F, G, drag, min);
        double z0 = getHorizontalComponent(B, D, vz, drag, min);
        return new Vector(x0, y0, z0);
    }

    private static double getHorizontalComponent(double posDiff, double D, double v, double drag, double time) {
        return D * (posDiff + time * v) / (pow(drag, time + 1) - 1);
    }

    private static double getVerticalComponent(double D, double X, double F, double G, double drag, double time) {
        double at = pow(drag, time);
        return (X * (1 - at) + time * F + G) / (D * (1 - at * (D + 1)));
    }

    private static double getNewtonRaphsonApproximation(double C, double O, double U, double N, double T, double E, double R, double P, double L, double A, double Y, double S, double drag, double guess, double accuracy, int maxIterations) {
        double f;
        double f1;
        double counter = 0;
        do {
            f = f(C, O, U, N, T, E, R, P, L, A, Y, S, guess, drag);
            f1 = f1(C, O, U, N, T, E, R, P, L, A, Y, guess, drag);
            guess -= (f / f1);
            counter ++;
        } while (Math.abs(f) > accuracy && counter < maxIterations);
        return guess;
    }

    private static double f(double C, double O, double U, double N, double T, double E, double R, double P, double L, double A, double Y, double S, double guess, double drag) {
        return C * pow(drag, 4 * guess)
            + O * pow(drag, 3 * guess)
            + U * pow(drag, 2 * guess)
            + N * pow(drag, guess)
            + T * sq(guess)
            + E * guess
            + R * sq(guess) * pow(drag, 2 * guess)
            + P * guess * pow(drag, 2 * guess)
            + L * guess * pow(drag, 3 * guess)
            + A * sq(guess) * pow(drag, guess)
            + Y * guess * pow(drag, guess)
            + S;
    }

    private static double f1(double C, double O, double U, double N, double T, double E, double R, double P, double L, double A, double Y, double guess, double drag) {
        double dragLog = log(drag);
        return C * 4 * dragLog * pow(drag, 4 * guess)
            + O * 3 * dragLog * pow(drag, 3 * guess)
            + U * 2 * dragLog * pow(drag, 2 * guess)
            + N * dragLog * pow(drag, guess)
            + T * 2 * guess
            + E
            + R * (2 * guess * pow(drag, 2 * guess) + pow(guess, 2) * 2 * dragLog * pow(drag, 2 * guess))
            + P * (pow(drag, 2 * guess) + guess * 2 * dragLog * pow(drag, 2 * guess))
            + L * (pow(drag, 3 * guess) + guess * 3 * dragLog * pow(drag, 3 * guess))
            + A * (2 * guess * pow(drag, guess) + pow(guess, 2) * dragLog * pow(drag, guess))
            + Y * (pow(drag, guess) + guess * dragLog * pow(drag, guess));
    }

    private static double sq(double x) {
        return x * x;
    }
}