!----------------------------------
!
!  MODULE ScheffTest  
!
!----------------------------------
!
!
! Module for testing the envelope space charge routines in
!   SADScript.  Specifically for testing  the functions found 
!   in the packages
!
!      Scheff.n
!      Trace3dToSad.n
!      TwissUtility.n
!
! Currently the file is set up to simulation the J-PARC transport
!   line at 181 MeV.
!
!    Author    : Christopher K. Allen
!    Created   : November, 2005
!




!!========================================
!!
!!  Initialize SAD
!!
!!========================================


FFS;	! Begin SADScript



!
!  GLOBAL CONSTANTS
!

strBeamline = "L3BT01all";        ! beamline 
strFileOut = "ScheffTestOut.txt"  ! output file name



!
!  	Load Beamline
!

GetMAIN["~ckallen/J-Parc/linac/simdb-LI_L3BT01-nopmq0000.sad"];
L3BT01 = ExtractBeamLine["L3BT01all"];


!
!	Initialize SAD Environment
!

USE L3BT01;
TRPT;
INS;
CAL;

NOCOD;
RFSW;

$DisplayFunction = CanvasDrawer;



!
!	Define the Initial Beam Particle
!

!MASS   = 0.939294 GEV;
!CHARGE = -1;
!MOMENTUM = 0.610624 GEV;

!InitialOrbits =  {{ 0.0,0.0,0.0,0.0,0.0,0.0 }};





!!========================================
!!
!!  FUNCTIONS
!!
!!========================================


!---------------------------------------
!
!  FUNCTION  SaveMatrix
!
!---------------------------------------
!
!
!  Function saving an arbitrary matrix to persistent storage.
!
!  Parameters
!    strFile   file name to store matrix
!    mat       matrix to be stored
!
!  Returned Value
!    None
!
!
!    Author  :  Christopher K. Allen
!    Created :  November, 2005
!

SaveMatrix[strFile_, mat_] := Module[
    {
    dims,          ! matrix dimensions vector
    M,             ! number of rows
    N,             ! number of columns
    m,             ! loop control - rows
    n,             ! loop control - columns
    fos            ! file output stream 
    },

    dim = Dimensions[mat];
    M = dim[[1]];
    N = dim[[2]];

    fos = OpenWrite[strFile];
    Write[fos, "Matrix Dimensions ", M, "x", N];

    For[m=1, m<=M, m++,
        For[n=1, n<=N, n++, 
 
           WriteString[fos," : ",  mat[[m,n]]]
        ];
        Write[fos, " : "];
    ];
    Close[fos];
];






!
!  INITIALIZE SIMULATION
!


!
!    TRACE3D Parameters
!

f  = 324.0e6;        ! RF frequency (Hz)

Er = 939.29432e6;    ! particle rest energy (eV)
W  = 181.0338e6;     ! beam kinetic energy (eV)
XI = 30.0e-3         ! beam current (A)

vecTwissXt3d = {-0.44117, 5.774, 1.889};
vecTwissYt3d = {0.21808, 6.4229, 1.706};
vecTwissZt3d = {0.3095, 2.0888, 466.99};



!
!    Compute SAD Parameters
!

Q  = XI/f;                            ! beam bunch charge (C)


vecTwissX = TraceToSadTransTwiss[vecTwissXt3d];
vecTwissY = TraceToSadTransTwiss[vecTwissYt3d];
vecTwissZ = TraceToSadLongTwiss[f, Er, W, vecTwissZt3d];



!
!  RUN SIMULATION
!

! Compute generalized perveance and initial moment matrix

K0   = ComputePerveance[f, Er, W, Q];   	
sig0 = CorrelationMatrix6D[vecTwissX, vecTwissY, vecTwissZ];


! Run simulation
{lstPos, lstGamma, lstSig} = ScheffSimulate[K0, sig0, 0.01];

! Store results
SaveBeamMatrixData[strFileOut, lstPos, lstGamma, lstSig];

!  Look at the Results
PlotBeamBeta[lstPos, lstSig];

!{lstPos, lstGamma, lstTm} = GetBeamlineElementData[];
!TmRf = lstTm[[161]];
!posRf = lstPos[[161]];
!SaveMatrix["SadRfGapMatrix.txt",TmRf];




Exit[];