(*
  Usage: Please set the fowllowing variables from outside as
  command line arguments, like this:
    wolframscript -file calculate_modes.wl psi4_file r rfilename ifilename
  for example:
    wolframscript -file calculate_modes.wl psi4.txt 2.5 rmodes.txt imodes.txt
*)

(* Parsing arguments *)
argv = Rest@$ScriptCommandLine;
(* Print[argv] *)
Check[
  psi4file = argv[[1]];
  (* {l, m, r} = Interpreter["Number"]/@argv[[2;;4]]; *)
	r = Interpreter["Number"]@argv[[2]];
	rfilename = argv[[3]];
	ifilename = argv[[4]];
,
  Print["Error: Check your command line arguments!"];Abort[]
];

(* Set Assumptions *)
$Assumptions = Pi > theta > 0 && 2 Pi >= phi >= 0 &&\
                Element[l, Reals] && Element[m, Reals] && -l<=m<=l;

Y[s_, l_, m_, theta_, phi_] :=
    (-1)^m*
      (Sqrt[(((l + m)!*(l - m)!)*(2*l + 1))/((((l + s)!*(l - s)!)*4)*Pi)]*
        Sin[theta/2]^(2*l))*
        Sum[
          (((Binomial[l - s, r]*Binomial[l + s, r + s - m])*(-1)^(l - r - s))*
          E^((I*m)*phi))*Cot[theta/2]^(2*r + s - m),
          {r, 0, l - s}
        ]

(* Y[l_, m_, theta_, phi_] := ((-1)^(l + m) E^(I m phi)
  Binomial[-2 + l, -2 - m] Cot[theta/2]^(-2 -m)
  Sqrt[((1 + 2 l) (l - m)! (l + m)!)/((-2 + l)! (2 + l)!)]
  Hypergeometric2F1[-2 - l, -l - m, -1 - m,
          -Cot[theta/2]^2] Sin[theta/2]^(2 l))/(2 Sqrt[Pi]) *)

psi4[theta_, phi_] = ToExpression[Import[psi4file, "Text"] ]/.{t->0};

mode[l_, m_] := Module[{f},
  (* f = Simplify[Conjugate[Y[-2, l,m,theta,phi] ]*psi4[r,theta,phi]*Sin[theta] ]; *)
  f = Conjugate[Y[-2, l,m,theta,phi] ]*psi4[theta,phi]*Sin[theta];
  Quiet[NIntegrate[f,{theta, 0, Pi}, {phi, 0, 2 Pi}] ]
]

(* mode[0,0,2.2]//Print *)

(* result = ParallelTable[mode[l,m],{l,1,8},{m,-l,l}]//Flatten; *)
lmlist = Flatten[Table[{l, m}, {l, 1, 8}, {m, -l, l}], 1]
result = ParallelMap[mode@@#&, lmlist]
rmodes = Re/@result;
imodes = Im/@result;
Export[rfilename, rmodes, "Table"]
Export[ifilename, imodes, "Table"]

(* result = mode[l,m]
Print[Re@result, " ", Im@result] *)