"""Reproduce conditional photo diagnostics; never changes model geometry or poses.

Run with the CAD Python environment (numpy/scipy/Pillow/matplotlib). Original
HEIC metadata was reconciled into the adjacent input JSON; only relevant tags
are retained. The original files are not needed to reproduce the diagnostics.
"""
from pathlib import Path
import json,math,hashlib
import numpy as np
from scipy.optimize import least_squares
from scipy.spatial.transform import Rotation
from scipy.optimize import minimize_scalar
from PIL import Image
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import fit_owner_side_view as side_fitter
from front_wing_geometry import shoulder_y,side_x
from front_clip_layout import panel_z
CAD=Path(__file__).parent;ROOT=CAD.parent;PATH=CAD/'owner-photo-registration-audit.json'
data=json.loads(PATH.read_text());engine=data['engine_photo'];side=data['side_photo']
held=['midget-1500-measured-core.json','midget-1500-measured-core.step','headlamp-reference.json','headlamp-reference.step','bonnet-structure-reference.json','bonnet-structure-reference.step','front-power-support-reference.json','front-power-support-reference.step','owner-side-view-6329.json','front_wing_geometry.py','bonnet_surface.py']
hashes={n:hashlib.sha256((CAD/n).read_bytes()).hexdigest()for n in held}
expected=data.get('results',{}).get('unchanged_geometry_and_camera_sha256',{})
if expected and any(hashes.get(n)!=value for n,value in expected.items()):
 raise RuntimeError('Historical photo audit source geometry changed. Preserve this baseline and use front-upper-profile-validation.json for the current revision.')

def focal(size,capture):
 metadata=next(x for x in data['captured_metadata']if x['source_basename']==capture)
 equivalent=float(metadata['FocalLengthIn35mmFormat'].split()[0])
 return equivalent*math.hypot(*size)/math.hypot(36,24)

def lip_outline():
 result=[]
 for x,z,begin in [(42,167,0),(-42,167,90),(-42,-167,180),(42,-167,270)]:
  for degree in np.linspace(begin,begin+90,17):
   a=math.radians(degree);result.append([x+18*math.cos(a),88,z+18*math.sin(a)])
 return np.array(result)

def project(parameters,points,f,principal):
 q=np.asarray(points)@Rotation.from_rotvec(parameters[:3]).as_matrix().T+parameters[3:6]
 return q[:,:2]/q[:,2,None]*f+principal

def edge_distances(samples,poly):
 d=np.roll(poly,-1,axis=0)-poly;delta=samples[:,None,:]-poly[None,:,:]
 t=np.clip(np.sum(delta*d[None,:,:],axis=2)/np.sum(d*d,axis=1),0,1)
 return np.sqrt(np.min(np.sum((delta-t[:,:,None]*d[None,:,:])**2,axis=2),axis=1))

size=np.array(engine['working_size']);principal=size/2;f=focal(size,'IMG_6855.HEIC')
points=np.array(engine['native_reference_points_mm']);observed=np.array(engine['observed_points_pixels']);lip=np.array(engine['lower_lip_samples_pixels']);outline=lip_outline();weights=np.array(engine['residual_normalization_weights_pixels'])
fits={}
for family in ['mixed_reference','lip_only_upper_points_held_out']:
 def residual(p):
  projected=project(p,points,f,principal)
  anchors=((projected-observed)/weights[:,None]).ravel() if family=='mixed_reference' else (projected[-1]-observed[-1])/weights[-1]
  return np.r_[anchors,edge_distances(lip,project(p,outline,f,principal))/6]
 candidates=[]
 for height in [650.,1000.,1500.]:
  position=np.array([150.,height,900.]);forward=np.array([0,100,0])-position;forward/=np.linalg.norm(forward)
  right=np.cross(forward,[0,1,0]);right/=np.linalg.norm(right);down=np.cross(forward,right);r=np.stack([right,down,forward])
  initial=np.r_[Rotation.from_matrix(r).as_rotvec(),-r@position]
  fit=least_squares(residual,initial,max_nfev=2500)
  depth=outline@Rotation.from_rotvec(fit.x[:3]).as_matrix().T+fit.x[3:6]
  if fit.success and np.isfinite(fit.x).all()and np.min(depth[:,2])>0:candidates.append(fit)
 assert candidates,'No finite positive-depth diagnostic fit'
 best=min(candidates,key=lambda x:x.cost);r=Rotation.from_rotvec(best.x[:3]).as_matrix();errors=project(best.x,points,f,principal)-observed
 fits[family]=dict(parameters=best.x.tolist(),focal_pixels=f,camera_engine_reference_mm=(-r.T@best.x[3:6]).tolist(),upper_point_residual_vectors_pixels=errors[:3].tolist(),upper_point_radial_rms_pixels=float(np.sqrt(np.mean(np.sum(errors[:3]**2,axis=1)))),lip_distance_rms_pixels=float(np.sqrt(np.mean(edge_distances(lip,project(best.x,outline,f,principal))**2))))
fig,ax=plt.subplots(figsize=(11,8),dpi=150);photo=Image.open(ROOT.parent/engine['image'].removeprefix('/cars/')).resize(tuple(size));ax.imshow(photo)
for family,color,label in [('mixed_reference','#00d9ff','Current source anchors + lip'),('lip_only_upper_points_held_out','#ff63d9','Lip only; upper points held out')]:
 p=np.array(fits[family]['parameters']);q=project(p,outline,f,principal);ax.plot(*np.vstack([q,q[0]]).T,color=color,lw=1.3,label=label);upper=project(p,points[:3],f,principal);ax.scatter(*upper.T,color=color,marker='x',s=35)
ax.scatter(*lip.T,color='#ffe66d',s=10,label='Hand-picked lip samples');ax.scatter(*observed[:3].T,edgecolor='black',color='#ffe66d',s=25,label='Approximate nut/cap observations');ax.set_xlim(350,1000);ax.set_ylim(1030,270);ax.legend(fontsize=8);ax.set_title('IMG_6855 · fixed EXIF-derived focal hypothesis; geometry unchanged',fontsize=10);fig.tight_layout();fig.savefig(ROOT/'body-side-study/engine-reference-audit.png');plt.close(fig)
# Keep the authoritative side camera untouched. Compare an additional EXIF branch.
side_data=json.loads((CAD/side['camera_source']).read_text());body=json.loads((ROOT/'body-data.json').read_text());target_f=focal(side['image_size'],'IMG_6329.HEIC');original_solver=side_fitter.least_squares
try:
 def constrained(fun,x0,**kwargs):return least_squares(lambda p:np.r_[fun(p),(p[6]-math.log(target_f))*10000],x0,**kwargs)
 side_fitter.least_squares=constrained;exif_branch=side_fitter.fit(side_data,body)
finally:side_fitter.least_squares=original_solver
samples=np.array(side['far_upper_wing_outline_samples_pixels']);zs=np.linspace(.23,1.99,300)

def side_project(camera,drop=0.):
 pts=np.array([[side_x(shoulder_y(z)-drop,z),shoulder_y(z)-drop,panel_z(z)]for z in zs]);q=pts@np.array(camera['world_to_camera_rotation']).T+np.array(camera['world_to_camera_translation'])
 return q[:,:2]/q[:,2,None]*camera['focal_length_pixels']+np.array(side['image_size'])/2

side_results={};fig,ax=plt.subplots(figsize=(12,6),dpi=150);ax.imshow(Image.open(ROOT.parent/side['image'].removeprefix('/cars/')))
for name,camera,color in [('retained_wheel_camera',side_data['camera_fit'],'#00d9ff'),('exif_focal_sensitivity',exif_branch,'#ff63d9')]:
 q=side_project(camera);errors=np.interp(samples[:,0],q[:,0],q[:,1])-samples[:,1]
 def cost(drop):
  p=side_project(camera,drop);return np.mean((np.interp(samples[:,0],p[:,0],p[:,1])-samples[:,1])**2)
 probe=minimize_scalar(cost,bounds=(0,.20),method='bounded')
 side_results[name]=dict(focal_pixels=camera['focal_length_pixels'],wheel_center_rms_pixels=float(np.sqrt(np.mean(np.array(camera['center_residuals_pixels'])**2))),upper_outline_y_residuals_pixels=errors.tolist(),upper_outline_y_rms_pixels=float(np.sqrt(np.mean(errors**2))),uniform_drop_probe_mm=float(probe.x*2032/2.4),probe_residual_rms_pixels=float(math.sqrt(probe.fun)),probe_applied=False,camera=camera)
 ax.plot(*q.T,color=color,lw=1.8,label=name.replace('_',' '))
ax.scatter(*samples.T,color='#ffe66d',edgecolor='black',s=30,label='Owner far upper-body outline');ax.set_xlim(1000,1710);ax.set_ylim(650,430);ax.legend(fontsize=8);ax.set_title('IMG_6329 · same front-wing geometry under two conditional cameras',fontsize=11);fig.tight_layout();fig.savefig(ROOT/'body-side-study/front-profile-audit.png');plt.close(fig)
for name,old in hashes.items():assert hashlib.sha256((CAD/name).read_bytes()).hexdigest()==old,name
# These diagnose reconstruction ambiguity; no residual threshold certifies fit.
data['results']={'engine_camera_families':fits,'front_profile_camera_sensitivity':side_results,'unchanged_geometry_and_camera_sha256':hashes,'engine_or_body_pose_applied':False,'native_geometry_changed':False,'interpretation':'The source-cover reference choices do not supply a calibrated common engine/body frame. The far front-body outline remains above owner samples under both side-camera hypotheses, identifying upper-profile reconstruction as a substantial next priority. The uniform-drop probes are diagnostics only: their spread and remaining residuals do not define fitted heights or justify moving the whole body/powertrain.','next_modeling_priority':'Reconstruct the longitudinal front upper-body profile and coordinate bonnet, wings, scuttle/hinges and closing interfaces. Preserve the nominal axles, marker and rear-body references; use existing owner views to test candidate geometry.','remaining_unknowns':['Fitted cover dimensions, hardware datum heights/stations and exact point correspondence','Principal point, digital crop, residual distortion and effective focal calibration','Common surveyed engine/body/mount axes and safe operating/removal clearance','Near-side reflection stripe is not a verified seam; far outline correspondence remains conditional']}
PATH.write_text(json.dumps(data,indent=2)+'\n')
print(json.dumps({'engine_diagnostics':{k:{a:b for a,b in v.items()if a in ['upper_point_radial_rms_pixels','lip_distance_rms_pixels']}for k,v in fits.items()},'body_diagnostics':{k:{a:b for a,b in v.items()if a in ['wheel_center_rms_pixels','upper_outline_y_rms_pixels','uniform_drop_probe_mm','probe_residual_rms_pixels']}for k,v in side_results.items()},'unchanged_geometry_and_camera_sources':len(hashes),'geometry_changed':False},indent=2))
