compiletime load "module_ppp.so";
case = BATCH.file;

PPP obj = {
  # Create forward reference for the maximal 
  # absorption at left and right side from center line:
  float AMAX : 0.5*(A_left_max + A_right_max + abs(A_left_max - A_right_max));
};

T0 = 15;     # Start time for time averaging
T1 = 20;     # End time for time averaging
DT = 2.5;    # Plot interval
radius = 75; # Target radius
tilt   = 10; # Magnet tilt angle

# Compute resulting center coordinate of absorption profile
pi = 4*atan(1);
Y_center = 92.5 - radius*sin(tilt*pi/180);

# Make center coordinate also visible in obj 
obj.Y_center = Y_center;

# Read averaged absorption and cell surfaces
for(time=T0; time<=T1+1e-6; time=time+2.5) {
  obj.read("$case/absorption/absorption_Arplus_$(time%,2f)us.pos", "A");
}
obj.read("$case/cell_surfaces.pos", "S_All");

# Consider only surfaces, where Arplus absorption is non-zero
obj.create_scalar("S", S_All*A>0 ? S_All ! 0);

# Create product A*S
obj.create_scalar("AS_product", S*A);

# Make horizontal cutplane at Target T2 averaging across upper half surface
obj.plane(-298, 0, 205,   298, 0, 205,    -298, 170, 205,   40, 0, 0, "average");
obj.reduce();

# Create averaged absorption sum A*S / sum S (if sum S > 0)
obj.create_scalar("A_avg", S>0 ? AS_product/S ! 0);

# Save plane view with projected absorption for illustration
obj.write_pos_scalar_binary("extract/$case-Arplus2D-T2.pos", 
                            "Averaged Arplus absorption %", 
                            A_avg/A_avg_max*100);
                          
# Make cumulated line profiles

# Create scalars for distinguishing absorption on left and right side from center
obj.create_scalar("A_left",  y<Y_center ? A_avg ! 0);
obj.create_scalar("A_right", y>Y_center ? A_avg ! 0);

# Convert plane profile into accumulated linescan
obj.xz_plane(Y_center, radius, "sum");

# Write relative line profiles into TXT file
obj.write_txt("extract/$case-T2-lineprofiles.txt", 
             ["x", "A_left_%", "A_right_%"], 
             [x, 100*A_left/AMAX, 100*A_right/AMAX]);