#!/bin/bash
until [[ `echo x$1` == 'x' ]]
do
    case $1 in
	-h) shift;
	    echo 'A tool to convert a particle data file to PeTar style, which can be used to start a simulation.'
	    echo 'Usage: petar.init [options] [particle data filename]';
	    echo '   particle data file: A file containing information of particles.';
	    echo '      The file contains 7 columns, representing mass, 3D positions and 3D velocities, respectively.';
	    echo '      Each line contains the data of one particle.';
	    echo '      If there are N binaries, they should be located in the first 2*N lines with two components next to each other.'
	    echo 'Options (default arguments shown in parentheses at the end):';
	    echo '  -f [S] Specify the output file (PeTar input data) name (default: input file name + ".input")';
	    echo '  -i [I] Skip the given number of rows in the input data file (default: 0)';
	    echo '  -s [S] Add stellar evolution columns: base | bse | no (default: no)';
	    echo '  -m [F] Set the mass scaling factor from the input data unit to [Msun]: mass[input unit]*m_scale=mass[Msun] (default: 1.0)';
	    echo '      Note that Msun is used as the mass unit in BSE based stellar evolution.';
	    echo '  -r [F] Set the radius scaling factor from the input data unit to [pc] (default: 1.0)';
	    echo '  -v [F] Set the velocity scaling factor from the input data unit to [pc/myr]  (default: 1.0)';
	    echo '          If "kms2pcmyr" is given, convert velocity unit from [km/s] to [pc/myr].';
	    echo '  -u     Calculate the velocity scaling factor based on the mass and radius scaling, then convert the data unit to (Msun, pc, pc/myr).';
	    echo '          The input data should use the Henon unit (total mass=1, G=1).';
	    echo '          The mass and radius scaling can be modified using -m and -r, respectively.';
	    echo '  -t     Add an external potential column and the position and velocity offsets to all particles in the header line.';
	    echo '         This is required when the external potential (e.g., Galpy) is enabled (--with-external in configure).';
	    echo '  -c [S] Set values of position and velocity offsets [in input unit], values are separated by "," (default: 0,0,0,0,0,0).';
	    echo '         The units will be transformed based on the scaling options (-r, -v, -u).';
	    echo '         This is required when the external potential (e.g., Galpy) is enabled and the option "-t" is used';
	    echo '  -R [S] Set the initial stellar radius for "-s base" mode (default: 0.0)';
	    echo 'Important notes:'
	    echo '    1) When using stellar evolution (e.g., BSE), be cautious with the scaling factor.';
	    echo '       It is recommended to use the unit set [Msun, pc, pc/myr] for the input data.';
	    echo '       If velocities are in km/s, use "-v kms2pcmyr" to convert the unit to pc/myr.';
	    echo '    2) If the input data is in the Henon unit, specify "-m" and "-r"; the option "-u" can calculate the correct velocity scaling factor.';
	    echo '    3) A PeTar-style input file generated by this tool can also be used as a data file to change scaling factors or stellar evolution columns.';
	    echo '       In this case, use "-i 1" to remove the header line.';
	    echo '       Data in BINARY format cannot be read.';
	    echo '    4) Be aware that this tool does not work for a PeTar snapshot from the middle of a simulation where binary positions may have changed.';
	    echo '    5) If an SSE/BSE-based stellar evolution package is compiled (e.g., --with-interrupt=bse), use the option "-s" to set the correct format of columns.';
	    echo '    6) If an external potential (e.g., Galpy) is used, use the option "-t" to set the correct format of columns.';
	    echo '       Additionally, use "-c" to set the position and velocity offsets of the cluster center in the galactic coordinate system.';
	    exit;;
	-f) shift; fout=$1; shift;;
	-i) shift; igline=$1; shift;;
	-s) shift; seflag=$1; shift;;
	-m) shift; mscale=$1; convert=1; shift;;
	-r) shift; rscale=$1; convert=1; shift;;
	-v) shift; vscale=$1; convert=1; shift;;
	-u) henon_unit=1; shift;;
	-R) shift; radius=$1; shift;;
	-c) shift; cm=$1; shift;;
	-t) extflag='yes'; shift;;
	*) fname=$1;shift;;
    esac
done

if [ ! -e $fname ] | [ -z $fname ] ; then
    echo 'Error, file name not provided' 
    exit
fi
[ -z $fout ] && fout=$fname.input
[ -z $igline ] && igline=0
[ -z $seflag ] && seflag='no'
[ -z $extflag ] && extflag='no'
[ -z $rscale ] && rscale=1.0
[ -z $mscale ] && mscale=1.0
[ -z $vscale ] && vscale=1.0
[ -z $radius ] && radius=0.0
[ -z $convert ] && convert=0
[ -z $henon_unit ] && henon_unit=0
[ -z $cm ] && cm='none'

echo 'Transfer "'$fname'" to PeTar input data file "'$fout'"'
echo 'Skip rows: '$igline
echo 'Add stellar evolution columns: '$seflag

n=`wc -l $fname|awk '{print $1}'`
n=`expr $n - $igline`

cm_array=''
if [[ $cm != 'none' ]]; then
    cm_array=(`echo $cm|sed 's/,/ /g'`)
    echo 'cm offset: pos: '${cm_array[0]}' '${cm_array[1]}' '${cm_array[2]}' vel: '${cm_array[3]}' '${cm_array[4]}' '${cm_array[5]}
#    awk -v ig=$igline -v x=${cm_array[0]} -v y=${cm_array[1]} -v z=${cm_array[2]} -v vx=${cm_array[3]} -v vy=${cm_array[4]} -v vz=${cm_array[5]} '{OFMT="%.15g"; if(NR>ig) print $1,$2+x,$3+y,$4+z,$5+vx,$6+vy,$7+vz}' $fname > $fname.off__
#else
#    awk -v ig=$igline '{if(NR>ig) print $1,$2,$3,$4,$5,$6,$7}' $fname >$fname.off__
fi

# first, scale data
if [ $convert == 1 ] && [ $henon_unit == 1 ] ; then
    G=0.00449830997959438
    echo 'Convert Henon unit to Astronomical unit: distance scale: '$rscale';  mass scale: '$mscale';  velocity scale: sqrt(G*ms/rs);  G='$G
    awk -v ig=$igline -v rs=$rscale -v G=$G -v ms=$mscale 'BEGIN{vs=sqrt(G*ms/rs)} {OFMT="%.15g"; if (NR>ig) print $1*ms,$2*rs,$3*rs,$4*rs,$5*vs,$6*vs,$7*vs}' $fname>$fout.scale__
    if [[ x$cm_array != x ]]; then
	cm_array_scale=(`echo ${cm_array[@]} | awk -v rs=$rscale -v G=$G -v ms=$mscale 'BEGIN{vs=sqrt(G*ms/rs)} {OFMT="%.15g"; print $1*rs,$2*rs,$3*rs,$4*vs,$5*vs,$6*vs}'`)
	cm_array=(${cm_array_scale[@]})
	echo 'cm offset (scaled): pos: '${cm_array[0]}' '${cm_array[1]}' '${cm_array[2]}' vel: '${cm_array[3]}' '${cm_array[4]}' '${cm_array[5]}
    fi
    mscale=1.0 # use for scaling from Petar unit to stellar evolution unit, since now two units are same, set mscale to 1.0
elif [ $convert == 1 ]; then
    [ $vscale == 'kms2pcmyr' ] && vscale=1.022712165045695
    echo 'Unit convert: distance scale: '$rscale';  mass scale: '$mscale';  velocity scale: '$vscale
    awk -v ig=$igline -v rs=$rscale -v vs=$vscale -v ms=$mscale '{OFMT="%.15g"; if (NR>ig) print $1*ms,$2*rs,$3*rs,$4*rs,$5*vs,$6*vs,$7*vs}' $fname >$fout.scale__
    mscale=1.0 # use for scaling from Petar unit to stellar evolution unit, since now two units are same, set mscale to 1.0
    if [[ x$cm_array != x ]]; then
	cm_array_scale=(`echo ${cm_array[@]} | awk -v rs=$rscale -v vs=$vscale '{OFMT="%.15g"; print $1*rs,$2*rs,$3*rs,$4*vs,$5*vs,$6*vs}'`)
	cm_array=(${cm_array_scale[@]})
	echo 'cm offset (scaled): pos: '${cm_array[0]}' '${cm_array[1]}' '${cm_array[2]}' vel: '${cm_array[3]}' '${cm_array[4]}' '${cm_array[5]}
    fi
else
    awk -v ig=$igline '{if(NR>ig) print $LINE}' $fname >$fout.scale__
fi
#rm -f $fname.off__

#         m,  r,        v,        bdata
base_col='$1, $2,$3,$4, $5,$6,$7, 0,' 
#         rs, id,    mbk, stat, rin, rout, acc_s, pot_t, pot_s, 
soft_col='0,  NR-ig, 0,   0,    0,   0,    0,0,0, 0,     0,'
if [[ $extflag == 'yes' ]]; then
    echo 'Add the external potential column (pot_ext)'
    soft_col=$soft_col' 0,' # pot_ext
    # header
    if [[ x$cm_array != x ]]; then
	echo '0 '$n' 0 '${cm_array[@]} >$fout
    else
	echo '0 '$n' 0 0.0 0.0 0.0 0.0 0.0 0.0' >$fout
    fi
else
    #header
    echo '0 '$n' 0' >$fout
fi
soft_col=$soft_col' 0' # nb

if [[ $seflag != 'no' ]]; then
    #       radius,  dm, t_record, t_interrupt
    se_col=$radius', 0,  0,        0,' 

    if [[ $seflag == 'base' ]]; then
	echo "Interrupt mode: base'
	echo 'Stellar radius (0): " $radius
	awk '{OFMT="%.15g"; print '"$base_col$se_col$soft_col"'}' $fout.scale__ >>$fout
    elif [[ "$seflag" == *"bse"* ]]; then
	#       type, m0,  m,     rad, mc,  rc,  spin, epoch, time, lum
	bse_col='1, $1*ms, $1*ms, 0.0, 0.0, 0.0, 0.0,  0.0,   0.0,  0.0,'

	echo 'Interrupt mode: '$seflag
	echo 'mass scale from PeTar unit (PT) to Msun (m[Msun] = m[PT]*mscale): ' $mscale
	awk -v ms=$mscale  '{OFMT="%.15g"; print '"$base_col$se_col$bse_col$soft_col"'}' $fout.scale__ >>$fout
    else
	echo 'Error: unknown option for stellar evolution: '$seflag
    fi
else
    awk '{OFMT="%.15g"; print '"$base_col$soft_col"'}' $fout.scale__ >>$fout
fi
rm -f $fout.scale__
