#Purpose: # To become familiar with VPython and prepare for what I will spend the summer of 2013 # working on. # I have been told that I will have to write a code in python that simulates a photon # entering a roman pot detector due to the Cherenkov radiation of a proton. The photon # is to enter a rectangular peice of quartz, reflect off a 45 degree boundary, and # then travel up another rectangular peice of quartz. # # Author Date Status # ================ ================ ============================ # Alexander Sage June 06, 2013 Original code: incomplete from visual import * scene.height = 800 #set size of window in which to run scene.width = 800 #the simulation #define the photon as a sphere with following characteristics photon = sphere(pos=(4,-30, 2), radius=1, color=color.yellow, material=materials.emissive) photon.vel = vector(2, 5, 2) temp = vector(0, 1, 0) #assign placeholder variable # set the coordinates of numerous solid rectangles to be used as # 3D walls backwall = box(pos=(-25, -30, 0), size=(0.5, 10, 10), color=color.white) #walls flat in x-y plane bottomside1 = box(pos=(0, -30, -5), size=(50, 10, 0.5), color=color.cyan) bottomside2 = box(pos=(0, -30, 5), size=(50, 10, 0.5), color=color.cyan, opacity=0.2) #walls flat in x-z plane bottomside3 = box(pos=(5, -35, 0), size=(60, 0.5, 10), color=color.cyan) bottomside4 = box(pos=(0, -25, 0), size=(50, 0.5, 10), color=color.cyan) #45 degree boundary between rectangles slope = Polygon([(-25, -35), (-35, -35), (-35, -25)]) extrusion(pos=[(0, -0, -5), (0, 0, 5)], shape=slope, color=color.red) #walls flat in x-y plane topside1 = box(pos=(30, -5, -5), size=(10, 60, 0.5), color=color.blue) topside2 = box(pos=(30, -5, 5), size=(10, 60, 0.5), color=color.blue, opacity=0.2) #walls flat in y-z plane topside3 = box(pos=(35, -5, 0), size=(0.5, 60, 10), color=color.blue) topside4 = box(pos=(25, 0, 0), size=(0.5, 50, 10), color=color.blue) topwall = box(pos=(30, 25, 0), size=(10, 0.5, 10), color=color.white) t = 0 dt = 0.005 while t < 100: #starting a while loop to update the position of the photon #and restrict its position to inside the quartz rate(1000) #how fast time travels photon.pos = photon.pos + photon.vel*dt #updating position # the sphere has a radius of 1 :: these if statements are written so that the photon # will reflect once the edge of the photon 'touches' one of the boundaries. #horizontal reflection if photon.pos.x - 1 <= backwall.pos.x: photon.vel.x = -photon.vel.x if photon.pos.x + 1 >= topside3.pos.x: photon.vel.x = -photon.vel.x if ((photon.pos.x - 1 <= topside4.pos.x) and (photon.pos.y >= -25)): photon.vel.x = -photon.vel.x #vertical reflection if ((photon.pos.x <= 25) and (photon.pos.y + 1 >= bottomside4.pos.y)): photon.vel.y = -photon.vel.y if photon.pos.y - 1 <= bottomside3.pos.y: photon.vel.y = -photon.vel.y if photon.pos.y + 1 >= topwall.pos.y: photon.vel.y = -photon.vel.y #depth reflection if photon.pos.z - 1 <= bottomside1.pos.z: photon.vel.z = -photon.vel.z if photon.pos.z + 1 >= 5: photon.vel.z = -photon.vel.z #slope reflection if photon.pos.y <= photon.pos.x + pow(2, 0.5) - 60: temp.x = photon.vel.x photon.vel.x = photon.vel.y #slope reflection swaps the x and y photon.vel.y = temp.x #components of velocity t = t + dt #Execution: #this code has an error when the photon hits the line (25, -25, 5) -> (25, -25, -5) #it travels inside the top wall of the horizontal rectangle until it exits through #the line it entered from. I don't have any ideas left how to fix it. #set initial conditions for photon.pos=(4, -30, 2) photon.vel = vector(8, 5, 2) #let run for t > 30 and see for yourself #it works fine for reflection from all surfaces, no other vertices were tested.