| 123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263 |
- import json, math, hashlib
- from pathlib import Path
- import numpy as np
- from PIL import Image, ImageFilter
- BASE=Path('.scratch/cockpit-shell-panel-layout/verify/artifacts/t35')
- PILOT=BASE/'negative-estate-pilot/input-v2'
- OUT=BASE/'roof-imagery-quality'
- OUT.mkdir(parents=True,exist_ok=True)
- configs=json.loads((PILOT/'roof-textures.json').read_text())
- regions=json.loads((BASE/'roof-quality-pilot/report.json').read_text())
- def pixels(path):
- return np.asarray(Image.open(path).convert('RGB'),dtype=np.float32)
- def stats(arr):
- rgb=arr/255.0
- lum=.2126*arr[:,:,0]+.7152*arr[:,:,1]+.0722*arr[:,:,2]
- sat=arr.max(2)-arr.min(2)
- dx=np.abs(np.diff(lum,axis=1));dy=np.abs(np.diff(lum,axis=0))
- return {'sizePx':[int(arr.shape[1]),int(arr.shape[0]),'source channels RGB'],
- 'luminancePercentiles':{str(p):round(float(np.percentile(lum,p)),2) for p in [1,5,25,50,75,95,99]},
- 'darkPixelFractionLumaUnder55':round(float(np.mean(lum<55)),5),
- 'brightPixelFractionLumaOver210':round(float(np.mean(lum>210)),5),
- 'saturationPercentiles':{str(p):round(float(np.percentile(sat,p)),2) for p in [50,90,99]},
- 'edgeGradientP95':round(float(np.percentile(np.concatenate([dx.ravel(),dy.ravel()]),95)),2),
- 'rgbStddev':np.round(np.std(arr.reshape(-1,3),axis=0),2).tolist()}
- textures=[]
- tile_map={}
- for bid,meta in configs.items():
- arr=pixels(meta['path']);rec={'id':bid,'roofCrop':stats(arr),'sha256':hashlib.sha256(Path(meta['path']).read_bytes()).hexdigest(),
- 'sourceTileCount':len(meta['sourceTiles']),'sourceTiles':meta['sourceTiles']}
- textures.append(rec)
- for t in meta['sourceTiles']:tile_map[t['path']]=t['sha256']
- for target in regions['results']:
- for src in target['processed']['sourceTiles']:tile_map[src['path']]=src['sha256']
- # Tile-edge continuity proxy. Compare RGB values immediately beside adjacent source tile boundaries.
- tiles={p:pixels(p) for p in tile_map}
- seams=[]
- items=sorted(tiles)
- for p in items:
- parts=Path(p).parts
- try:z=int(parts[-3]);x=int(parts[-2]);tms=int(Path(p).stem)
- except (ValueError,IndexError):continue
- y=2**z-1-tms
- for axis,nx,ny in [('east',x+1,y),('south',x,y+1)]:
- q=f'D:/GIS/GIS/yuanlingPic3/{z}/{nx}/{2**z-1-ny}.png'
- if q not in tiles or p>=q:continue
- a,b=tiles[p],tiles[q]
- edgeA=a[:,-1,:] if axis=='east' else a[-1,:,:]
- edgeB=b[:,0,:] if axis=='east' else b[0,:,:]
- d=np.sqrt(np.mean((edgeA-edgeB)**2,axis=1))
- la=.2126*edgeA[:,0]+.7152*edgeA[:,1]+.0722*edgeA[:,2]
- lb=.2126*edgeB[:,0]+.7152*edgeB[:,1]+.0722*edgeB[:,2]
- seams.append({'tileA':p,'tileB':q,'orientation':axis,'meanRgbRmsDelta':round(float(np.mean(d)),2),
- 'p95RgbRmsDelta':round(float(np.percentile(d,95)),2),'meanLumaDelta':round(float(np.mean(np.abs(la-lb))),2)})
- report={'generatedAt':'2026-10-15','method':'Pillow/Numpy raster statistics only; not roof segmentation, radiometric correction or facade evidence',
- 'roofCrops':textures,'sourceTileCount':len(tiles),'adjacentTileSeams':seams,
- 'interpretationLimits':['Dark-pixel share and luminance spread combine cast shadows, roof material, capture mosaics and building occlusion.','Adjacent tile edge mismatch is a proxy and may cross different land cover; it is not a seam detector.','The 512/205px GLB atlases are roof bounding rectangles; only the roof mesh UV subset is rendered.','No source image is altered; no geometric claim follows from these statistics.'],
- 'runtimeConfigChanged':False,'runtimeRendered':False}
- (OUT/'roof-imagery-statistics.json').write_text(json.dumps(report,indent=2))
- print(json.dumps({'passed':True,'roofCrops':[(x['id'],x['roofCrop']['darkPixelFractionLumaUnder55'],x['roofCrop']['luminancePercentiles']['50'],x['sourceTileCount']) for x in textures], 'sourceTileCount':len(tiles),'seams':[(x['orientation'],x['meanRgbRmsDelta']) for x in seams]}))
|
