Instructions 1. Submit two files for each problem set. The first is a R Markdown (.Rmd) file that can be run without error from start to end. The second is a PDF rendered from your R Markdown file or created using LATEX. 2. Name your files following this convention: [Last Name]_ps1.Rmd and [Last Name]_ps1.pdf. 4. Comment on your code wherever possible and explain your ideas in detail. You will get credit for showing the steps you take and for explaining your reasoning, even if you do not get the correct final result.Part 1: Fixed Effects and Random Effects Models In this exercise, we will use the dataset sibling_data.dta to study the effect of mother experiencing stress in pregnancy on a child’s birth weight. All variables are described in the table below. Variable Name Variable Detail Group Variable householdid Unique household id; Siblings from the same household share the same household id Dependent Variable birthwt birth weight measured in pounds Independent Variable stress Whether or not mother experienced stress in pregnancy (Yes=1; No=0) female Infant’s gender (Female=1; Male=0) magebirth Mother’s Age at Birth numsibling Number of Sibling meduy Mother’s Years of Schooling feduy Father’s Years of Schooling 1.[20pts] Import the dataset sibling_data.dta to your R environment. Including all independent variables, build an OLS model, a fixed effects model, and a random effects model (you can use the plm R package). Summarize regression results in a table. Make sure to label your models in your table (the column.labels argument in stargazer()). 1 Hint: householdid is used as a grouping variable in your FE or RE model. 2.[5pts] Interpret the coefficient of stress in the OLS model. 3.[5pts] Interpret the coefficient of stress in the fixed effects model. 4.[5pts] Interpret the coefficient of magebirth in the fixed effects model. 5.[5pts] Perform a F-test to compare the OLS and the fixed effects models and interpret outputs (use pFtest from the plm R package). 6.[5pts] Compared to an OLS model, what are the benefits of using a fixed effects model? 7.[5pts] Why don’t we get coefficients of meduy, feduy, and numsibling in the fixed effects model? 8.[10pts] Why do we get coefficients of meduy, feduy, and numsibling in the random effects (random intercepts) model? 9.[10pts] How does the coefficient of stress change between (a) the OLS model and the fixed effects model and (b) the OLS model and the random effects model? What could be the potential causes that lead to these changes? 10.[10pts] Perform a Hausman test to compare fixed-effects and random-effects models (figure out how to use the phtest function from the plm R package; you might look into the help page and examples by typing ?plm::phtest into your console) (a) What is the null hypothesis of this model? What is the alternative hypothesis? (b) Is the null hypothesis rejected? Based on the test, which model would you use? Hint: Hausman test evaluates whether the results (i.e. the estimated coefficients) from a fixed-effects and random-effects model are significantly different. You are encouraged to search online what the test is and how to decide which model is preferred. 2
Contents Instructions 1 Prerequisite 1 Part 1 Multivariate Regression 2 Part 2 Logistic Regression 2 Part 3 Causality 3 Instructions 1. Submit two files for each problem set. The first is a R Markdown (.Rmd) file that can be run without error from start to end. The second is a PDF rendered from your R Markdown file or created using LATEX. 2. Name your files following this convention: [Last Name]_ps3.Rmd and [Last Name]_ps3.pdf (for example, Jiang_ps3.Rmd). 4. You are encouraged to discuss the problems with your classmates. But the R Markdown and PDF files that you submit have to be created on your own. Please do not ask for solutions from students in earlier cohorts. 5. Comment on your code wherever possible and explain your ideas in detail. You will get credits for showing the steps you take and for explaining your reasoning, even if you do not get the correct final result.Prerequisite Load multiple packages to your environment. knitr::opts_chunk$set(echo = TRUE) ## load packages here library(dplyr) ## ## Attaching package: ‘dplyr’ ## The following objects are masked from ‘package:stats’: ## ## filter, lag ## The following objects are masked from ‘package:base’: ## ## intersect, setdiff, setequal, union Part 1 Multivariate Regression gss.csv contains a dataset of the following variables: Variable Name Variable Detail id A respondent’s unique ID abortion A respondent’s answer to the question “Should a Woman be Able to Have Legal Abortion”; -1 = No, 0 = It depends, and 1 = Yes pid A respondent’s party identification; 0 = Strong Democrat, 1 = Democrat, 2 = Democrat Leaning, 3 = Independent, 4 = Republican Leaning, 5 = Republican, 6 = Strong Republican linc A respondent’s income (logged) female A dummy variable of sex; 1 = Female, 0 = Male mar A dummy variable of marital status; 1 = Married, 0 = Not married 1. [5pts]Plot the distribution of abortion using a bar plot • Hint: treat both abortion and pid as continuous variables and fit regular OLS 3. [5pts]Calculate heteroskedasticity-robust standard errors for Model 2; call this model Model 3. 4. [5pts]Create a regression table with three columns: Model 1, Model 2, and Model 3. 5. [10pts]Plot the predicted value of abortion of male and female respondents with different pid values. Be clear about which variables you fix at what values when generating the plot. 6. [5pts]Using both the tables and the plot, interpret the results. Part 2 Logistic Regression 1. [5pts]Replicate the linear probability model and logistic model fitted in Lab 8, then fit a new logistic regression model using the same set of predictors but treat pid as an ordered categorical variable. • Hint: by ordered categorical variable, I mean you need to factorize pid, or create dummy categories 2. [5pts]Create a regression table with three columns: Model 1 (linear probability model), Model 2 (logistic regression with pid as a numeric variable), and Model 3 (logistic regression with pid as an ordered categorical variable). Part 3 Causality A study on COVID-19 constructed a “COVID risk factor” score based on the COVID infection rate of a given area (defined by zip code). A researcher wants to estimate the effect of having a vaccination center in the area on that area’s COVID risk factor score. She compiled a dataset that contains each area’s COVID risk factor score and whether the area has a vaccination center. She then estimated the effect of having a vaccination center using the “naive estimator” we discussed in class. You noted that the quality of information residents have about COVID and the vaccine can be a confounding variable that affects both the area’s infection rate and whether there is a vaccination center in the area. Assume that you are able to estimate the relationships this “informedness” confounder (info) and the original “vaccination center” predictor (vaccine) have with the COVID risk factor score (covid_risk), which can be simulated using the following code (n is sample size): library(“scales”) set.seed(1234) # set the same seed to ensure identical results e = rnorm(n, 0, 0.5) covid_risk = rescale( 0 – 7*vaccine – 2*info + e, to = c(0, 100)) 1. [5pts]Import the data covid.csv. According to the potential outcome framework, create two new variables y_c and y_t where y_c is the value of the potential outcome of “risk factor” when the individual is not treated, and y_t is the value of the potential outcome of “risk factor” when the individual is treated. Explain your steps. 2. [5pts]Fill out the table below (round to 1 decimal points): Group YT YC Treatment Group (D = 1) E[YT|D = 1] =? E[YC|D = 1] =? Control Group (D = 0) E[YT|D = 0] =? E[YC|D = 0] =? 3. [15pts]Estimate/calculate the following (write down the formula/equation you use): i. The Naive Estimator of ATE ii. The “true” Average Treatment Effect iii. The “true” Treatment Effect on the Treated iv. The “true” Treatment Effect on the Control v. Selection Bias
Instructions 1. Submit two files for each problem set. The first is a R Markdown (.Rmd) file that can be run without error from start to end. The second is a PDF rendered from your R Markdown file or created using LATEX. 2. Name your files following this convention: [Last Name]_ps2.Rmd and [Last Name]_ps2.pdf (for example, Jiang_ps2.Rmd). 4. You are encouraged to discuss the problems with your classmates. But the R Markdown and PDF files that you submit have to be created on your own. Please do not ask for solutions from students in earlier cohorts. 5. Comment on your code wherever possible and explain your ideas in detail. You will get credits for showing the steps you take and for explaining your reasoning, even if you do not get the correct final result.knitr::opts_chunk$set(echo = TRUE) ## load packages here Part 1: Assumptions Consider the “true” data-generating process (DGP): Y = β0 + β1X + 1. [10pts] Suppose we obtain a sample of size n from this process. State the assumptions under which the OLS estimator βˆ1 is BLUE (Best Linear Unbiased Estimator) for β1. 2. [10pts] For each of the assumptions, discuss what will go wrong when the assumption is violated. 3. [5pts] Let β0 = −0.25, β1 = 1.2, X ∼ Γ(5,4), and ∼ Normal(0,1). Here, Γ(α,ψ) denotes the Gamma distribution with shape parameter α and rate parameter ψ. A simple Google search will lead you to a way to simulate from this distribution (Don’t get scared away from the new distribution and the fancy Greek letters!). Simulate a dataset of size n = 3,000 from this process in which all of the assumptions you’ve discussed above hold. 4. [20pts] For each of the assumptions, create a plot which illustrates how the violation of the assumption affects the regression results (all other assumptions should be satisfied). This can be a scatter plot, a plot of the regression line, the sampling distribution of the OLS estimator (comparing what the regression table suggest with actual simulations), or anything else. The important point is to show how the violation leads us to false decisions if we assume the assumption is true, not just showing that the assumption is violated. When simulating data, you don’t have to use the parameters set in the previous problem. Also, you can skip the “no autocorrelated errors” assumption, which goes beyond the course materials. Further bigger hint: For the violation of exogeneity assumption, the easiest way is to create two variables, x1 and x2, with non-zero covariances between the two. y is a function of x1 and x2, while when you specify the OLS regression, you omit either x1 or x2. In this case, the error term, by construction, will be correlated with the independent variable. Then discuss the difference between the true model, where exogeneity assumption is satisfied. Part 2 Simulate and Plot Correlation [15pts] Replicate the following figure which demonstrates the different levels of Pearson’s correlation. Hint: Look up mvrnorm() from the MASS package.Part 3 Simulate Interactions 1) [15pts] Simulate a data generation process that will provide you data to roughly replicate the two-panel figure below. You should show your code and replicate the two-panel example figure. Note: Your output does not need to replicate the exact layout of the example figures. You will get credit as long as you generate a similar plot that illustrates the relationship clearly. Remember to use set.seed() for any random process.2) [15pts] Simulate a scatter plot that demonstrates the Simpson’s Paradox:3) [10pts] Choose one of the above data simulations you just generated, write down the regression equation you estimated using OLS, and interpret the coefficients in your own words. Hint: You will need to specify group-specific slopes and intercepts. Part 4: The Replication Project Read through the article and determine which variables you will need to replicate the project. Go to the IPUMS webpage and submit a screenshot of the page where you can see the samples and variables that you have selected (you can obtain this from your “Data Cart”).
Contents Instructions 1 Prerequisite 1 Part 1. Functions 2 Part 2. Population, sample, and sampling distribution 2 Part 3. Simulate confidence interval 3 Instructions 1. Submit two files for each problem set. The first is a R Markdown (.Rmd) file that can be run without error from start to end. The second is a PDF rendered from your R Markdown file or created using LATEX. 2. Name your files following this convention: [Last Name]_ps1.Rmd and [Last Name]_ps1.pdf (for example, Jiang_ps1.Rmd). 4. You are encouraged to discuss the problems with your classmates. But the R Markdown and PDF files that you submit have to be created on your own. Please do not ask for solutions from students in earlier cohorts. 5. Comment on your code wherever possible and explain your ideas in detail. You will get credits for showing the steps you take and for explaining your reasoning, even if you do not get the correct final result.Prerequisite Load multiple packages to your environment using the following code (you can add more packages to the current list as per your need): knitr::opts_chunk$set(echo = TRUE, fig.pos = ‘H’, fig.align=”center”) library(tidyverse) library(kableExtra) library(gridExtra) Part 1. Functions Recall the formulas for population mean: and variance: 1 N µ = Xxi N i=1 (1) 1 N σ2 = X(xi − µ)2 (2) N i=1 where N is the population size. 1.[5pts] Write a function that calculates the population mean according to Equation 1 without using any R functions that directly calculate the mean. For example, you cannot use mean() from base R, or summarise(., mean = mean()) from tidyverse. 2.[5pts] Write a function that calculates the population variance according to Equation 2 without using any R functions that directly calculate the variance. For example, you cannot use var() from base R, or summarise(., var = var()) from tidyverse. You cannot use mean() from base R, either. 3.[10pts] Import gapminder.csv to your R environment. • Apply the two functions you just created to the lifeExp variable in gapminder. • Use base R functions that directly calculate mean and variance to the same lifeExp variable vector. • Report your results of the above two steps either in text or in a table. The results for the mean should be equal, but the results for variance should be different. Find out and explain why the results in variance differ. Part 2. Population, sample, and sampling distribution 1.[3pts] Create a population that follows a normal distribution with population mean µ = 5 and population variance σ2 = 1 with 100,000 observations. 2.[3pts] Create a histogram of the population with appropriate title and labels. Add a vertical line at the population mean. 3.Draw a random sample from the population, with sample size n = 50. 4.[3pts] Plot a histogram of the sample with appropriate title and labels. Add a vertical line at your point estimate of the population mean. How does this histogram compare to the one you created in question 2? 5.[6pts] Based on your sample, report your point estimate, µˆ, of the population mean and the standard error of this estimate. 6.[6pts] Simulate the sampling distribution of the sample mean (n = 50) using 1000 draws. That is, repeat the action you took for question 3 for 1000 times and save the mean you get for each repetition to a numeric vector. 7.[3pts] Create a histogram of the sampling distribution of the sample mean you simulated in question 6 with appropriate title and labels. Add a vertical line at your point estimate of the population mean. 8.[6pts] Using the sampling distribution you obtained in question 6, report your point estimate of the population mean µˆ and the standard error of this estimate. 9.[20pts] Repeat questions 3 to 8 increasing the size of your sample to n = 1000. Plot and report your results. Then, using the concepts that we learned in class, summarize the differences with respect to what you obtained with a sample of 50 (Hint: there are two related concepts, one related to the sample estimate of the population mean, and the other related to the uncertainty of the standard error of this estimate). Part 3. Simulate confidence interval 1.[10pts] In your own words, explain what is a 95% confidenece interval and what is a p-value. 2.[5pts] Write down the mathematical formulas you use to calculate the 95% confidence interval of the sample mean for a given sample (assume the sample size is larger than 50 so you can use the z instead of the t-statistics). Make sure you explain the meanings of the notations you use in your formulas. 3.[15pts] Replicate the plot below and interpret the plot in one or two sentences. Hint: The plot shows the simulation result when a sample (size n) is randomly drawn 100 times without replacement from a hypothetical population with mean µ and standard deviation σ. µ has a correct, fixed answer, but n and σ depend on your choice, as long as the simulation results look similar to the plot. Note: As there is randomness involved in the simulation process, your plot does not need to replicate the exact same result as the example. Also, try your best to tweak the layout of your plot, but it does not need to look exactly the same as the example. You will get credits as long as you generate a similar enough graph that convey the conceptual points clearly and correctly. Remember to use set.seed() for any random process.
Logistics • Quiz explanations Part 1: Replication Project Tips 1.1 Samples • For the 1990 sample, use the 1% metro sample • For the 2010 sample, use the single-year ACS sample, not 3- or 5-year pooled sample • If you read pp.1046 carefully, you will notice that Rs with the top and bottom earning percentile are excluded – You can create percentiles using quantile(ma$WEEKEARN, seq(0.01,1,0.01)), suppose your dataframe is ma, and the weekly earning variable is WEEKEARN 1.2 Variables • Use BPL rather than NATIVITY • The latter has no valid values for most samples • Use HISPAN to exclude Hispanic Whites and Hispanic Blacks • Use CLASSWKR to determine whether R is in a public sector or not – You should look at CLASSWKRD, which gives detailed classification of CLASSWKR • Use CPI99 to adjust inflation for INCWAGE • The main dependent variable is the logged form of weekly earnings – You will need WKSWORK1 and WKSWORK2 to measure the number of weeks worked last year. WKSWORK1 always gives the best continuous estimate, but when WKSWORK1 is not available, you should turn to WKSWORK2 • To estimate potential years of experience, the formula is given by LMEXP = AGE – EDUYEAR – 6 – EDUYEAR needs to be estimated – Codes for this process are available in the code folder 1.3 Duncan’s Dissimilarity Index • In Table A1a and A1b, you will notice that there is a dissimilarity index. This is a very commonly used measure of occupational segregation. – The standard Duncan’s Dissimilarity/Segregation Index is given by: 1 X D = n ai − bi 2 A Bi=1 where ai and bi is the number of White and Black workers in occupation group i. A and B represents the total number of White and Black workers. • Instead of using hundreds of OCC1990 categories, you will use 2-digit aggregated categories of OCC1990 – Codes are available on Brightspace – Use the merge() function 1.4 General Instructions • It is totally okay if you cannot get exactly the same numbers! I also couldn’t. • But they should be close enough. If they deviate a lot, you need to explain your speculations why the numbers differ this much. Part 2: Causality: The Potential Outcome Framework 2.1 The Fundamental Problem of Causal Inference • The modern way of thinking about causality is to think about outcomes in a counterfactual approach • For example, the effect of a policy treatment on an outcome Y , is to think about the difference between Yit, i.e., the potential outcome of individual i receiving the treatment, and Yic, i.e., the potential outcome of the same individual, if not receiving the treatment. Either one of the two terms is never observed. – Notation-wise, Y c and Y t are both potential outcomes (i.e., Rubin’s approach) – ATE is defined as ATE = E[Yit] − E[Yic] = E[Yit − Yic] – We can only observe E[Yit|Di = 1] and E[Yic|Di = 0] 2.2 Naive Estimation of the Average Treatment Effect • At the population level, the average treatment effect (ATE) is defined as: τ = ATE = E[Yit − Yic] = E[Yit] − E[Yic] * Since we do not observe the population level Y T or Y C, the naive approach to estimate the population level ATE uses the following equation: τˆ = E[Yit|Di = 1] − E[Yic|Di = 0] • which calculates the difference in the expected value of Yi in the observed treated group (E[Yit|Di = 1]) and the expected value of Yi in the observed control group (E[Yic|Di = 0]). • The estimated naive ATE will be unbiased if the assignment to treatment is purely random. 2.3 Selection Bias • However, if there are selection bias that lead to certain kinds of unit to go into the treatment or control group, the naive estimator will be biased. τˆ = E[Yit|Di = 1] − E[Yic|Di = 0] = E[Yit|Di = 1] − E[Yic|Di = 1]+E[Yic|Di = 1] − E[Yic|Di = 0]| {z } | {z } ATT selection bias where E[Yit|Di = 1]−E[Yic|Di = 1] is the treatment effect on the treated and E[Yic|Di = 1]−E[Yic|Di = 0] is the selection bias. You can think of it as the baseline difference of Yi if both the treatment and the control group are not treated. • There is also a definition of the treatment effect on the control (ATC), which can be expressed as E[Yit|Di = 0] − E[Yic|Di = 0] • E[Y ic|Di = 1,Xi] − E[Yic|Di = 0,Xi] = 0| {z } selection bias – This is a strong assumption, i.e., strong ignorability assumption – Failure to satisfy the assumption will lead to omitted variable bias, or selection bias, or the violation of zero-conditional mean assumption Part 2 Exercise Assuming you know both potential outcomes Yit and Yic on the same individual, as well as their realized outcomes. Answering the following questions:1. Calculate ATE: 2. Calculate ATT: 3. Calculate ATC: 4. Naive estimate of the ATE: 5. What causes the naive ATE to deviate from the true ATE in this example? Part 3: Some Examples of Causal Inferential Studies • The most straightforward approach is to manipulate treatment assignment to be completely random. In this case, selection bias is zero, by design. • Effect of neighborhood stigma on economic transactions, The Proceedings of the National Academy of Sciences (PNAS) – Random assignment of seller’s neighborhood information on an online market • The mark of a criminal record, AJS – Random assignment of racial background of job seekers with criminal record • Exposure to Opposing Views can Increase Political Polarization: Evidence from a Large-Scale Field Experiment on Social Media, PNAS – Random assignment of Republican or Democrats voters to follow twitter accounts from elected officials and opinion leaders with opposing political views • In many cases, experiments with random assignment is either unethical or unfeasible. People therefore turn to natural experiments, where assignment of treatment is quasi-random. Nation building through foreign intervention: Evidence from discontinuities in military strategies, QJE • Quasi-random assignment of bombing and airstrikes in the Vietnam War Lifetime earnings and the Vietnam era draft lottery: evidence from social security administrative records, AER • Quasi-random assignment of military participation in the Vietnam War • Instrumental variables are seen as quasi-experiments if their variations come from certain social or natural designs • Commitment through Sacrifice: How Longer Ramadan Fasting Strengthens Religiosity and Political Islam, ASR – IV from the time-shifting feature of Ramadan that makes the fasting duration vary from year to year • Community and the Crime Decline: The Causal Effect of Local Nonprofits on Violent Crime, ASR – IV from nonprofit organizations in art, media, and medical industries to instrument local nonprofit organizations targeting violent crimes
library(tidyverse) library(stargazer) library(coefplot) library(sandwich) library(lmtest) library(ggpubr) Part 0: Housekeeping • There are multiple ways to include a figure in Rmarkdown. Please see one example below. knitr::include_graphics(“graph/multivariate_types.png”)Figure 1: Types of Multivariate Relationships Part 1: How Does Multivariate Relationships Affect Regression Estimates • Multiple Causes ## empty results woz
Part 1: Multivariate Regression & Interaction with One Dummy Dummies • For categorical variables, we create dummies or convert them to 0 or 1 dummies when we want to include them in a regression model • Note that for a categorical variable that have n categories, the regression model will only have n − 1 dummies or categorical variable predictors, because the nth dummy is redundant given that if an observation does not belong to any of the n − 1 category, then it must be in the nth category • We call the left-out category the reference category • Question: what if we include all n categories? • You should always interpret your model coefficients with the reference category in mind. This could get complicated when you have multiple dummy variables, especially when they are interacted in your model In the case of the dummies representing “race” in the earnings_df that we will be using today, we have: Category Dummy1(black) Dummy2(other) White 0 0 Black 1 0 Other 0 1 Exercise 1. Import earnings_df.csv to your environment. Perform the following data cleaning steps: (1) If age takes the value 9999, recode it as NA; (2) Create a new variable female that equals 1 when sex takes the value female, and equals to 0 otherwise; (3) Create a new variable black that equals 1 when race is black and equals to 0 otherwise; (4) Create a new variable other that equals to 1 when race is ’other‘ and 0 otherwise. 2. Use the describe() function from the psych package to generate a quick descriptive statistics of your data. 3. Now, estimate the following models and display your model results in a single table using stargazer(m_1, m_2, …, m_n, type=”text”). (1) Model 1: earn ~ age (baseline) (2) Model 2: earn ~ age + edu (3) Model 3: earn ~ age + edu + female (4) Model 4: earn ~ age + edu + female + race (5) Model 5: earn ~ age + edu + female + race + edu*female 4. Write down your prediction equation for Model 5 Ii = β0Ai + β1Ei + β2Fi + β3Bi + β4Oi + β5Ei × Fi 5. In Model 5, holding other variables constant, what will be the predicted difference in estimated mean earnings for a white man and a white women? 6. Holding other variables constant, what will be the predicted difference in estimated mean earnings for a white women and a black women? White woman = β0A¯i + β1E¯iWhite woman = β0A¯i + β1E¯i + β2 + β5E¯i 7. Holding other variables constant, what will be the predicted difference in estimated mean earnings for a white man and a black woman? White man = β0A¯i + β1E¯iBlack woman = β0A¯i + β1E¯i + β2 + β3 + β5E¯i earnings_df
SummaryFor this lab you will implement a basic LC3 assembly program which calculates the sum of all odd numbers stored between x3200 and x32FF inclusive.Program DetailsThe final sum will be stored at memory location x3300. Once stored, your program will halt the machine using the HALT command.You must place any .FILL pseudo-ops after your HALT command. Placing .FILL pseudoops before the HALT will cause the data stored at those locations to be executed as instructions. This can cause your code to produce unexpected results and prevent Gradescope from executing correctly.Bonus: +10 points for making and using an IsEven subroutine to determine if a specific value is even or odd. The subroutine must start at a label named IsEven and expect a single input parameter (a 2s complement number) in R0. If the given number is even, the subroutine will return 0. If the number is odd, it will return any non-zero value. The result of the subroutine will be stored in R0. You must define your subroutine after your halt instruction so that the subroutine will only execute if jumped to.The IsEven subroutine must also implement callee save and preserve the values of registers R1-R5. NOTE: If your subroutine doesn’t touch a register (i.e. modify its value) then you don’t need to push that registers value onto the stack!Development NotesWhen testing if a number is odd or even, please remember that a number can be negative or positive and to account for that in your program. Also note that if you use the user stack, you will need to initialize R6 with the stack pointer value first. The user stack should always start at xFDFF per the LC3 ISA.For developing and testing your assembly program I recommend using the online LC3 simulator at: https://wchargin.com/lc3web/ . While this is the simulator I recommend, you are allowed to use any LC3 simulator to develop this project. **I also highly recommend downloading the LC3’s ISA from Canvas to use as a reference for this project. It is an invaluable resource which is vital for the project.** When saving your assembly program for Gradescope, you will need to use a text editor such as sublime text or notepad++AND R0,R0,#0A good comment would be: ;sets R0 to 0A poor comment would be: ;ANDs the value of R0 to 0 and stores the resultFinal NotesRubricCategory Points Autograder 70 Code Formatting 15 Comments 15
Fall Semester Telecommunications Circuits Laboratory EPFL6 Fading Channels In the previous Section 5 we have seen that the signal undergoes a random phase-rotation in the channel. In practical radio communication systems there is an additional effect a communication systems engineer has to take into account. The magnitude of the signal varies randomly over time. This effect is called fading. The reason for fading can be divided in two groups: Large-scale fading is already addressed in our model by means of the SNR value. Throughout this exercise we are only interested in dealing with the small-scale fading case. You have probably already observed small-scale fading as large variation of the signal strength of your radio or cell phone. The presence of small-scale fading changes the input-output relation from the pure additive white Gaussian-noise (AWGN) case y[k] = x[k] + n[k] (6.1) to y[k] = h[k]x[k] + n[k] (6.2) where h is a (complex) random variable with certain statistical properties. You will notice that in small-scale fading environment the BER decays much slower with increasing SNR compared to the AWGN case. (Note: Despite the misleading naming there is also an additive white Gaussian noise component in the fading case). The AWGN channel can also be modeled as a fading channel with a deterministic h = 1. 6.1 Coherence The propagation conditions do not change instantaneously; hence, channel coefficients are usually correlated in time and frequency. This means that if the channel is in a bad condition, then neighboring (in time or frequency) channel coefficients are also very likely in a bad condition.This correlation of the channel coefficients is usually described by two parameters: Thecoherence time and the coherence bandwidth. These parameters describe the distance in time/frequency over which the channel changes significantly. While these measures look very intuitive there is no strict definition of how large or likely the change of the channel has to be. There are multiple different definitions of the term. So consider these numbers not as precise values, but as orders of magnitude. From the viewpoint of the receiver fading can be divided in two scenarios: • Slow fading, the change of the channel condition within your time-span of interest (e.g. one frame) is negligible, but nevertheless has to be taken into account between successive transmission or different channels. • Fast fading, the change of the channel within your time-span has significant influence on the performance of your system and has to be compensated in your receiver. 6.2 Diversity psuccess ≈ 1 − pout If the same data is sent over L independent channel realizations, then all coefficients have to be in outage for none of the copies to be received with proper signal strength. The probability of successful transmission is therefor much higher and becomes psuccess ≈ 1 − (pout)L 6.3 Receive Diversity It is possible to obtain diversity on the receiving side without any modifications of the transmitter. One can employ multiple antennas in order to receive the same signal multiple times under different fading conditions. The antennas have to be sufficiently spaced apart so that their distance exceeds the coherence distance. Therefore this is also called space diversity. There are two common schemes to make use of space diversity. Antenna Selection Diversity Antenna selection diversity is a cheap way to implement receive diversity. It is achieved by attaching two or more sufficiently spaced antennas with a switch to the same receiver. The receiver periodically compares the signal strength at the antennas and selects the strongest signal for reception. hequivalent = max(h1,h2,…) (6.3) This approach offers the advantage that it is very cheap in terms of hardware costs as only one RX chain is required. On the downside this scheme can not exploit the full possible diversity gain. Maximum Ratio Combining y = hx + n (6.4) A matched filter of the form h∗ hmf = ||h|| leads to an equivalent channel of: (6.5) h∗h h∗ hmf · y = x + n = ||h||x + n˜ ||h|| ||h|| (6.6) hequivalent = ||h|| (6.7) 6.4 Rayleigh Fading A channel model which is often used to describe the effect of fading on the channel is the Rayleigh fading. The channel coefficient is modeled as a complex circular-symmetric Gaussian random variable h ∼ CN(0,1). (6.8) 6.5 Your Tasks For the next exercise we assume a slow fading channel with Rayleigh Fading with three receive antennas. We assume the exaggerated case where the channel coefficient is constant in magnitude and phase during a frame, but completely independent from the coefficient of any other frame or antenna. To keep the slow fading assumption valid we can only transmit shorter frames. A6T1 Extend the frame synchronization code to also provide an estimate of the magnitude of the channel. (Note: Don’t forget that there is already some kind of normalization.) A6T2 Implement three receiver configuration and compare their BER: (1) Single Antenna: Use only the signal received from the first antenna (2) Antenna Select: Use the signal received with the maximum amplitude (3) Maximum Ratio Combining: Combine the signal receiver from all available antennas
Fall Semester Telecommunications Circuits Laboratory EPFL4 Timing Synchronization The timing problem is easily understood: The receiver samples the incoming signal with the known sampling rate 1/Ts, but how does it know exactly when to sample? It does not know this, it can only sample the signal with a random time offset, and it is the task of the time synchronization unit to estimate this time offset and correct it. 4.1 Effect of a Timing Offset The effect of the timing offset can be modelled by including the filter δ(t−εT) into the effective channel. The parameter ε denotes the time offset, normalized by the symbol duration T. The matched filter output then becomes . (4.1) Take a look at Figure 4.1, where several RC pulses are plotted. Imagine that we want to detect the data symbol a[0]. If we sample the signal exactly at t = 0, all pulses except of the middle one are zero, and the received sample is equal to the data symbol plus noise. Now assume that we have a timing error of, say, ε = 0.1, which means that the signal arrives a bit too late, and hence our receiver samples the signal at t = −0.1T. At this time, however, the center pulse has not yet reached its maximum amplitude; furthermore, the contributions of the other pulses do not vanish, which leads to inter-symbol interference. In order to attain ISI-free reception, we must attempt to generate symbols (4.2) where εˆ is an estimate of the timing error. If εˆ = ε, then z(nT + εTˆ ) = a[n] + w0(nT + εTˆ ) and the data can be optimally recovered. 4.2 Synchronization via Interpolation In principle, there are two possibilities to generate the desired symbols z(nT + εTˆ ). The first possibility is to adjust the time instant of A/D conversion. However, this results in an undesired feedback to the analog part of the system. We will therefore use a different method that worksFigure 4.1: Raised cosine pulses with α = 0.22Figure 4.2: Time synchronization unit 4.6. 4.3 The Estimation Algorithm The remaining question is how to obtain the estimate εˆ. For the derivation of the estimator we use the sequence x[k] = x(kTs) which is defined as the instantaneous signal power x(kTs) = |z(kTs)|2. (4.3) Appendix 4.5 shows the following relation: if there is no time offset (ε = 0), then the expectation of x(kTs) is a real-valued and even signal. Accordingly, its spectrum(i.e., the expectation of X(f)) is then real-valued for all frequencies f. Now, if we introduce a nonzero offset εT in z(kTs) and, accordingly, in x(kTs), the original real spectrum is multiplied by e−2πjεTf. Evaluating the spectrum at the symbol rate f = 1/T then results in a complex value whose phase is −2πε. We can thus use arg{X(1/T)} (4.4) as an estimator of the timing offset (see Appendix 4.5 for details). The spectral component X(1/T) can be calculated by performing a discrete Fourier transform (DFT). In general, the spectrum X(f) is given as , (4.5) where we have assumed that the signal starts at t = 0. The spectral component at f = 1/T is thus (4.6) arg . (4.7) The cumulative sum acts as an averager: During the reception, more samples become available, and the estimate becomes more and more reliable. It is easily seen that the algorithm (4.7) only works if the true parameter ε is constant during the whole reception, which is the case in this lab. In practice, however, the timing offset will change over time, and the estimation accuracy would be severely degraded if we took all available samples (even old ones that were caused by a different ε than the current) into account. An important measure here is the coherence time, which is the duration over which a parameter can assumed to be approximately constant. Let the coherence time for the timing offset be denoted as Tc;ε = Nc;εT. Then the estimation algorithm can be implemented as a sliding sumover the 4Nc;ε most recently received samples: . (4.8) 4.4 Your Tasks A4T1 Implement the estimation of the timing offset εˆ as described in (4.7). Note that for each new sample n, the sum will only include 4 new terms. These terms can be accumulated to the sum obtained for the sample n − 1. A4T2 Implement both a linear and a cubic interpolator as presented in the appendix 4.6. The value of εˆ·(oversampling factor), provides two indications: the integer part indicates between which samples to perform the interpolation, and the fractional parts indicates the value that needs to be interpolated. A4T3 Take advantage of the preamble samples to obtain a better estimate of εˆ. 4.5 Appendix In this section, we show that 4.4 is an asymptotically unbiased estimator of the timing offset ε.(An estimator is asymptotically unbiased if it becomes unbiased as the number of observations goes to infinity.) 4.5.1 Signal Model We assume M-PSK signaling, which means that the amplitude of the data symbols a[n] is constant. Without loss of generality, set |a[n]| = 1. The symbols are assumed to be uncorrelated,i.e. E{a[m]a∗[n]} = δ(m − n). The received signal after matched filtering is given as , (4.9) where gRC(t)is the combined impulse response of the pulse shaping filter and the matched filter,and is therefore real-valued and even. In our lab, gRC(t) is a raised cosine pulse. Furthermore, w0(t) is the AWGN with variance σw2 after matched filtering. For the digital realization, the signal is sampled at four times the symbol rate: T/Ts = 4. 4.5.2 Continuous Time Consideration We start with examining the noise-free signal z˜(t) = z(t) − w0(t). Let x¯(t) = E{|z˜(t)|2} (4.10) denote the expected value of the instantaneous signal power. Since the data symbols are uncorrelated, x¯(t) becomes (4.11) The function x¯(t) is periodic with period T with n0 = n − 1 (4.12) and can therefore be explained into Fourier series . (4.13) The Fourier coefficients are given as (4.14) The first Fourier coefficient, i.e. the spectral component X¯(1/T), is substitute u = t − εT (4.15) Thus, we can calculate the timing offset perfectly from the noise-free and time-continuous received signal: . (4.16) 4.5.3 Discrete Time Realization Now we consider the noisy and sampled signal x[k] = |z(kTs)|2 and show that (4.17) arg{Xˆ(1/T)} (4.18) is an asymptotically unbiased estimate of ε. We start by calculating the expected value of x[k]: E{x[k]} = E{|z˜(kTs) + w0(kTs)|2} = E{|z˜(kTs)|2}+E{z˜∗(kTs)w0(kTs)}+E{z˜(kTs)w0∗(kTs)}+E{|w0(kTs)|2}| {z } | {z } | {z } | {z } x¯(kTs) 0 0 σw2 = ¯x(kTs) + σw2 (4.19) with the noise-free signal x¯(t) that we have examined in the previous section. Assuming that the signal starts at k = 0, the corresponding spectrum is . (4.20) Evaluating this for f = 1/T yields ∞ E{X(1/T)} = Xx¯(kTs)e−2πjk/4. (4.21) k=0 If we approximateE{X(1/T)} at time t = nT by using the 4Nc;ε most recently received samples (the Nc;ε most recently received symbols), we get the estimate . (4.22) Since x¯(t) is bandlimited to ±2/T we can calculate this from x¯(t), rather than from the discrete x¯(kTs): (4.23) with c1 calculated in (4.15). We can now calculate the expected value of εˆ:argE{Xˆ(1/T)} (4.24)4.6 Appendix 4.6.1 Perfect Interpolator Let f(x) be a (possibly complex-valued) function of one real-valued parameter. Assume that we only know f(x)for x ∈ Z, i.e. a series of samples. If f is bandlimited with a maximum frequency of 1/2, then the function is completely described by the samples. The problem of finding f(x) for x 6∈ Z is called interpolation. The ideal interpolator is a lowpass filter, i.e. the function can be calculated from its samples as . (4.25) However, this interpolator is not realizable because it has infinitely many taps. 4.6.2 Linear Interpolator The simplest possible interpolator is the linear one. Assume without loss of generality that we are interested in the intermediate value f(x) with x ∈ (0;1). The linear interpolator then calculates f(x) using the two closest samples f(0) and f(1) according to f(x) = (1 − x)f(0) + xf(1). (4.26) The linear interpolator can be implemented very efficiently, but its accuracy can be problematic. Consider the situation illustrated in Figure 4.3. Given are the functionand √ the samples f(0) = f(1) = 1/ 2. Now, for all x ∈ (0;1), the linear interpolator returns √ f(x) = 1/ 2, which is obviously a very bad approximation of the true value if x is around 0.5. On the other hand, if we wanted to estimate a value with x ∈ (−1;0) using the samples f(−1) and f(0), the linear estimation would be quite good. We see that the linear interpolator yields especially bad results if the function has a local maximum or minimum between the two considered samples.Figure 4.3: Linear interpolator for 4.6.3 Cubic Interpolator The interpolation can be improved by taking one additional sample on each side into account, which leads us to the cubic interpolator. Again we assume that we want to estimate f(x) with x ∈ (0;1), this time using the samples f(−1),f(0),f(1) and f(2). The idea is to approximate f(x) with a cubic polynomial f˜(x) and then evaluate it at the desired point. In general, f˜(x) = ax3 + bx2 + cx + d, where the coefficients a,b,c,d must satisfy (4.27) . (4.28) With a matrix inversion we get the result . (4.29)
Fall Semester Telecommunications Circuits Laboratory EPFL5 Phase Synchronization In a previous lab we saw that the real part of the base band signal is modulated onto a cosine wave at radio frequency, and the imaginary part onto a sine wave. Since the cosine and sine functions are orthogonal, the two signal parts can be separated in the receiver. However, we run into a similar problem as the one discussed in the last lab: The receiver merely sees an oscillating waveform (a sinusoid): r(t) ∝ sin(2πfct + θ) during one symbol duration. Is it a sine, a cosine, or something in between? In order to differentiate between the two base functions, the receiver needs a phase reference, i.e. it needs to know the phase offset θ. The task of the phase synchronization unit is to provide this reference. 5.1 Signal Model The phase synchronization is performed after the interpolator and thus works on the timingcorrected symbols zε[n] = z(nT + εTˆ ). (5.1) We will assume a perfect timing recovery and interpolation, so that zε[n] = z(nT + εT). zε[n] = a[n]ejθ[n] + w0[n]. (5.2) A simple but realistic phase noise model is the random walk: The first phase offset θ[0] is uniformly distributed in [0;2π). Then, over time, the phase offset evolves according toFigure 5.1: Received symbols without (a) and with (b) phase offset θ[n] = θ[n − 1] + ∆θ[n] ( mod 2π) (5.3) where the difference between two successive phase offsets is a Gaussian random variable . (5.4) Figure 5.2 illustrates a typical phase noise process for a standard deviation of σ∆θ = 0.004 w 0.23◦. We see that the phase drift is rather slow, but still non negligible if the transmission is longer than about thousand symbols. Therefore, the phase offset must be tracked during the transmission. 5.2 Synchronization via Rotation The correction of the phase offset is straightforward. Since the phase offset causes a rotation of the signal space,zε[n] = a[n]ejθ[n] + w0[n]. , we simply rotate the signal space into the opposite direction: zε;θ[n] = zε[n]e−jθˆ[n] ˆ[n]) 00 (5.5) j(θ[n]−θ = a[n]e + w [n] with w00[n] = w0[n]e−jθˆ[n]. The phase synchronization unit is illustrated in Figure 5.3. 5.3 The Estimation Algorithm The estimation algorithm, which is known in the literature as the “Viterbi and Viterbi” algorithm, uses the fact that we have a QPSK modulation. We therefore know that (5.6)Figure 5.2: Typical phase noise process for σ∆θ = 0.004Figure 5.3: Phase synchronization unit for all transmitted symbols a[n], and thus . (5.7) The operation (·)4 removes the information from the data symbols. Exploiting this fact, we can derive the following estimator, which operates on a symbol-by-symbol basis (i.e. does not take any correlation between the phases into account): . (5.8) This estimator is approximately unbiased if the SNR is high: arg{−zε[n]4}} arg{−(a[n]ejθ[n] + w0[n])4}} (5.9) arg{ej4θ[n]} (c) = θ[n]. The approximation (a) is valid for a high SNR, because then the variation of the argument is small, and the nonlinear argument operation can be linearized. In this case, the expectation and the argument operation can be exchanged. Equation (b) holds because the noise has zero-mean and is uncorrelated with the data sequence. The step (c) does not allow us to uniquely identify any θ ∈ [0;2π). To see this, consider, for example, and (we drop the symbol index for simplicity). Then, the estimator will give us: (5.10) (5.11) We can only identify θ correctly if we know which quadrant it is in. One way to have access to this information is to consider the phase estimate for the previous symbol. Assuming that the phase changes slowly, the two consecutive phases will most likely be close to each other (BUT: not necessarily in the same quadrant). However, in order to use this method we need a initial phase information about which quadrant we can start from. 5.4 Initial Phase Estimation An initial phase estimate can be obtained by exploiting the preamble that is used for frame synchronization in the following way. Assume that there is no noise, so that z[n] = a[n]ejθ[n] , and let p[n] denote the preamble sequence of length Np . The frame synchronization function computes the following correlation: Np−1 c[n] = X p¯∗[i]z[n − i], (5.12) i=0 where p[i] = p[Np − 1 − i] is the reverted preamble sequence. c[n] is maximized when the sequences z[n − i] and p[i],i = 0,…,Np − 1, are equal. Let the n which maximizes the correlation be denoted by nmax . Assuming that θ[nmax − i] ≈ θ over the preamble, we have: Np−1 c[nmax] = X p¯∗[i]p[nmax − i]ejθ[nmax−i] (5.13) i=0 Np−1 ≈ X p¯∗[i]p[nmax − i]ejθ (5.14) i=0 Np−1 = ejθ X p¯∗[i]p[nmax − i]. (5.15) i=0 For the phase of c[nmax] we have: arg{c[nmax]} = θ. (5.16) Thus, the phase of the peak of the correlator can be used as an initial phase estimate. 5.5 Your Tasks A5T1 Create phase noise by implementing the random walk described in (5.3). A5T2 Modify the phase synchronization unit so that it provides an initial phase estimate according to (5.16). A5T3 Implement the phase synchronization unit and integrate it into your system. Compare your phase estimates with the true values. (1) Apply phase noise to signal with awgn noise. (2) Apply viterbi-viterbi algorithm on the payload data. (3) Correct the payload with the estimated phase error. (4) Consider a way to filter the phase according to the nature of random walk, write a new theta estimation trial in parallel with (2) and (3) to allow a comparison of the estimate. Note: MATLAB angle function returns angles in the interval [−π,π]. After step (3), the image will still experience some erroneous phase jump unless a very high SNR is used (e.g. 100dB).
Fall Semester Telecommunications Circuits Laboratory EPFL3 Pulse Shaping and Matched Filtering Up to now we have only considered the transmission of complex data symbols over an AWGN channel. Of course, in reality there is no such thing as a complex symbol coming out of the antenna. Instead, the symbols have to be transformed into a time-continuous waveform. This transformation is called modulation. 3.1 Modulation 3.1.1 Transmission at Radio Frequency (RF) The principle of digital communication systems is that only a finite number M of different symbols can be transmitted. Each symbol ai is mapped onto a pulse gi(t),i = 1,…,M. Since only finitely many pulses are possible, the receiver can recover the transmitted data if the signal is not too heavily distorted by the channel. In this lab course, we only consider linear modulation schemes, where the pulses gi(t) can be expressed as scaled versions of the same base pulse: gi(t) = aig(t). Hence, the magnitude of the data symbol determines the amplitude of the transmitted pulse, and the argument determines its phase. (An example of a nonlinear modulation scheme is frequency shift keying, where the information lies in the frequency of the transmitted signal.) Let T denote the duration of one data symbol, i.e., a new pulse is transmitted each T seconds. Note that the pulse g(t) can still be longer than T ; the pulses are allowed to overlap. The signal s(t) can then be written as the superposition of the pulses, scaled by the data symbols and delayed by multiples of T : . (3.1) The signal s(t) is called a complex baseband signal, because it is complex-valued and its spectrum is centered at f = 0. In order to transmit the signal over the air, it must be transformed into a real-valued signal and its spectrum must be shifted to a carrier frequency , where B is the bandwidth of the signal. This operation can be written as follows: sRF(t) =
Fall Semester Telecommunications Circuits Laboratory EPFL2 Frame Synchronization In the first lab, the received signal samples corresponded exactly to the transmitted symbols, i.e. the beginning and the end of the data was inherently given. In reality, we receive a signal continuously, and we have to detect the beginning of the data frames within the received signal. This operation is called Frame Synchronization. A common way to achieve frame synchronization is to prepend the data blocks with a special signal that is known to the receiver, called a preamble. The receiver can then search for the preamble using a correlation filter. If the preamble is long enough and has a random-like structure, it is very unlikely that the data signal (or the received noise) has a similar shape as the preamble, and therefore a peak in the correlator output indicates with high probability the beginning of a data frame. 2.1 The Need for a Random-Like Structure In the introduction above, we mentioned that the preamble should have a random-like structure. Why is this the case? Could we not instead use, e.g. a series of ones as the known data? After all, this would simplify the correlator, since all multiplications could be omitted and the correlator could simply be implemented as a summation over Np consecutive received symbols, where Np is the length of the preamble. First, it is likely that such a data sequence also occurs in the payload data. Consider for instance a black and white bitmap file, where a white pixel is stored as the symbol 0, and a black pixel as 1. Since there are many consecutive white or black pixels, the correlator would not be able to differentiate between the preamble and the actual data. By using a random-like data sequence as preamble, we can circumvent this problem.(a) (b) Figure 2.1: Autocorrelation of the all-ones sequence (a) and a random sequence (b) of length Np = 20.Figure 2.2: LFSR used for creating the frame synchronization sequence latter case, it will still be possible to identify the peak, as long as the SNR is not too low. It should be clear now that a random-like data sequence should be used as preamble. Why random-like? Obviously, the receiver needs to know the sequence in order to look for it, so the transmitter cannot use really random data. Hence, we need a sequence which has statistical properties similar to random noise, but which is still reproducible in the receiver. One possibility would be to create a random sequence once, and then store it in the transmitter and receiver. But this would be very memory consuming if the sequence is long. Another possibility is to use pseudo-noise (PN) sequences. These are sequences with randomlike statistical properties, but which can still be reproduced deterministically. A common way to generate PN sequences is by means of a linear feedback shift register (LFSR). The output of the LFSR is periodic and it can be shown that the statistical properties of the LFSR output closely resemble those of a really random sequence generated by flipping a fair coin Np times. In order to reproduce the sequence, the receiver only needs to the position of the LFSR feedback taps, the initial state of the LFSR, and the length of the sequence. In the lab, we use the LFSR depicted in Figure 2.2, initialized with all ones. The length of the preamble is Np = 100. The bits of the preamble are modulated using binary phase shift keying (BPSK), where the bit 0 is mapped to the symbol +1, and the bit 1 is mapped to the symbol −1. The Figure 2.3 illustrate the functionning of an LFSR over 4 iterations. n = 0 n = 1Figure 2.3: Illustration of Functionning principle of an LFSR 2.2 The Detection Algorithm Our goal is to detect the presence of the known preamble sequence p[i],i = 0,1,…,Np − 1 in the received signal r[n] = a[n] + w[n]. The main idea, as already mentioned, is to feed the received signal through a correlator, which correlates the received signal with the known preamble. As soon as the portion of the received signal that is inside the correlator closely resembles the preamble, the magnitude of the correlator output exhibits a distinct peak. How do we detect this peak? We could directly take the magnitude of the correlator output and compare it to a certain threshold value. However, the problem is that the absolute strength of the correlator output depends on the unknown SNR value. In order to use a fixed threshold for peak detection, the correlator output must be normalized with respect to the received signal power. The correlator output at discrete time n can be written as Np−1 c[n] = X p∗[i]r[n + i]. i=0 We decide on the presence of a peak if the following condition is fulfilled: (2.1) (2.2) The denominator in (2.2) is the energy of the received signal that is currently inside the correlation filter, and γ is the decision threshold. 2.3 Using the Profiler in MATLAB In the previous assignment, you used the commands tic and toc to measure the overall execution time of your code. For a more detailed analysis, you can use the MATLAB profiler, which gives you a detailed report including the execution of each line that was executed and how often it was called. The profiler can be accessed by using the following commands: • profile on starts the profiler, clearing previously recorded profiling statistics. No Transmission Preamble BPSK ImageData QPSK Figure 2.4: Signal with unknown number of time steps where no transmission occurs, followed by the preamble (BPSK) and the image data (QPSK). • profile off stops the profiler. • profile viewer stops the profiler and displays the results in the profiler window. • profsave saves the results in HTML format. The HTML files are stored in a subfolder of the current folder named profile results. 2.4 Your Tasks A2T1 Implement the LFSR shown in Figure 2.2 as a MATLAB function and generate the preamble sequence p[i]. A2T2 The matrix file task2.mat contains a signal where the start of the transmission occurs at an unknown time step and is structured as shown in Figure 2.4. (a) Implement the correlator as a MATLAB function following the equations (2.1-2.2). (b) Plot the normalized and unnormalized correlator output for several SNR values ranging from −5 dB to 10 dB. You should notice a distinct peak within the signal. A2T3 In this task, you will convert your correlator into a preamble detector. Write the detector function that returns the index of the first sample of the ImageData. (Hint: make use of the code of the correlator.) A2T4 In this task, you will read the data contained in the received signal, and display the picture transmitted. Your receiver should run the frame synchronization algorithm on the (noiseless) signal contained in task2.mat until the preamble is detected. Then, the transmitted image should be received and displayed. Since the signal is noiseless, if your choice of peak detection threshold γ is correct, the image should always be displayed perfectly. In order to highlight the importance of frame synchronization, pick a threshold that will result in false alarms (a very low threshold will do this). Bonus: Use the MATLAB profiler on your code. Which part of your code takes the longest time to run?
Fall Semester Telecommunications Circuits Laboratory EPFL1 Mapping and Demapping 1.1 The AWGN ChannelFigure 1.1: System model for Lab 1 We now replace the BSC by another very common but a bit more sophisticated channel model. In the additive white gaussian noise (AWGN) channel the symbols are complex numbers which represent the amplitude and phase of the transmitted signal and carry the digital data. After transmission over a channel which disturbs the symbols with additive noise, the demapper converts the received noisy symbols back to a bit stream. In the following, the transmitter and channel blocks are explained in more detail. • Source: Throughout the lab, we use grayscale images as payload data. The pixels are transmitted rowwise, starting with the upper left pixel. The brightness of each pixel is represented by an unsigned 8-bit value, and the least significant bit (LSB) is transmitted first. • Mapper: The modulation scheme used by our system is quadrature phase shift keying √ (QPSK), which means that the bits are pairwise mapped onto symbols a = (± ±j)/ 2, as shown in Figure 1.2. With this modulation scheme, the information lies solely in the phase, since all symbols have the same amplitude |a| = 1, hence the name Phase Shift Keying. • AWGN channel: In the channel, the transmitted symbols a[n] are disturbed by complexvalued additive white Gaussian noise (AWGN) w[n] with zero mean and variance . (A complex-valued Gaussian distributed random variable w is obtained by taking two independent real-valued Gaussian random variables with half the variance, w1,w2 ∼ , and combining them to a complex number w = w1 +jw2.) “White” in this context means that the noise samples are uncorrelated, i.e. E{w[m]w∗[n]} = σw2δ[m−n].Figure 1.2: Mapping of bit pairs b0b1 to QPSK symbolsFigure 1.3: BER for uncoded QPSK modulation, AWGN channel The received symbols are then given as r[n] = a[n] + w[n]. (1.1) An important measure in communication systems is the ratio of average signal power to average noise power (signal-to-noise ratio, SNR), usually expressed in decibel: SNR[dB] . (1.2)Figure 1.4: Alternative Mapping of bit pairs b0b1 to QPSK symbols Since the power of the data symbols is normalized to 1(|a[n]| ≡ 1), the SNR of our system is . The best achievable bit error rate (BER) for an uncoded QPSK modulation is plotted in Figure 1.3 as a function of the SNR. 1.2 Your Tasks In order to perform the labs on your own computer you need to install the Communications Toolbox. To install the toolbox on Matlab simply go to Home>Add-Ons>Get Add-Ons, search for Communications Toolbox and install it. A1T1 You are given a transmitted signal, as well as the width and height of the transmitted image. The signal can be loaded into the MATLAB workspace with the command load. Your task is to generate a noisy received signal by adding, dependent on a SNR value in dB, correctly scaled noise to the loaded signal.2 Use the MATLAB function randn to generate the noise. The received symbols are then demapped into bits and the received image is displayed. image decoder. Hint: Note that the signal length is not equal to the number of pixels in the image. opt. Plot the received symbols as clouds of points in the complex plane. Compare plots of the noisy constellations of the received signal for different channel conditions. A1T2 • For the signal space constellation diagram shown in Fig. 1.5, find the average energy of the red and blue constellations (independently).Figure 1.5: Constellation diagram • Normalize the average transmitted symbol energy to 1 for the two constellations inde- pendently. Plot the energy-normalized constellation points of the two schemes (red and blue) on the same figure. What observations can you make? A1T3 Now write a function that experimentally reproduces the BER plot from Figure 1.3. Create a random bit sequence, map it onto symbols according to Figure 1.2, and add scaled noise. Demap the noisy symbols to bits and determine the bit error rate. Repeat these steps for different SNR values and then create a BER plot with the function semilogy (in order to get a logarithmic scale on the y-axis). Finally, compute the BER plot where you use the symbol mapping as proposed in Figure 1.4. A1T4 An alternative image decoder based on lossy image compression is implemented by the function compressed decoder. Unfortunately the script was programmed by a lazy PhD student and is therefore full of errors. Your mission is to hunt down all bugs. Hints: • Keep track of the arrays size under the ”Workspace” subwindow in matlab. • Under the ”Run” button in matlab, you can choose to pause the code execution on errors, making the debugging easier. • Under the ”Run” button in matlab, you can set the code to execute when pressing the button, avoiding going back and forth between the main script and the function being debugged. • A breakpoint placed at the end of a function let you investigate the state of the variables just before they are freed by the function return. • When founding a mistake, make sure that it was the only location in the code where this error occurs. 1.2.1 Plotting in MATLAB Guidelines for plotting graphs: • Plot the BER in a logarithmic scale (semilogy). • Make sure that the ranges of the axis are meaningful. • Always label the axis (xlabel / ylabel). • Assure that each plot has a meaningful legend. (legend) • Use the MATLAB function for export of pictures (saveas) NOT a screenshot.
1 Backpropagation with logistic loss [40 points] Let fw : RD → R be a single-hidden-layer, fully connected, neural network with no biases and parameters w = {W(1),w(2)}, such that its forward pass computes x (1) . (2) Here, σ(·) denotes a sigmoid activation. In what follows, we will assume D = 5 and K = 6. 1. [5 points] Implement a function predict(X,W) that given a batch of B samples (an array X of size B×D) and a dictionary of weights (W={w 1: D x K, w 2: K x 1}) computes the set of outputs . Do not use for loops. 2. [5 points] Implement a function logistic loss(y, y hat) that given a vector of labels y ∈ {0,1}B and a vector of predictions yˆ ∈ RB computes the average logistic loss of a batch. Compute the average loss for the case yˆ = 0 and y = 0. Explain your result. Reminder : Given a true label y ∈ {0,1} and a soft prediction ˆy ∈ [0,1], the logistic loss is defined as L(y,yˆ) = −y log(ˆy) − (1 − y)log(1 − yˆ). (3) 3. [6 points] Implement another function stable logistic loss(y, z 2) that given a vector of labels y ∈ {0,1}B and a vector of logits z(2) ∈ RB computes the average logistic loss of a batch. Make sure that your function is stable. Compute the average loss for the case when z2 = −1010 · 1 and y = 0. Explain your result. Hint : You can use np.logaddexp. 4. [12 points] Derive analytically, using the backpropagation algorithm, the expressions for the partial derivatives of the loss with respect to the weights, i.e., and , (4) when the loss is computed using a single pair (x,y). Recall that we are using the logistic loss, and that the final implementation should be stable. 5. [12 points] Implement a function gradient(X,y,W) which computes the gradient (i.e., vector of partial derivatives) of the average loss with respect to all the weights for a batch . Do not use for loops. 2 Classifying FashionMNIST using neural networks [30 points] In this exercise you will play with the dataset called FashionMNIST. Some examples are shown in figure 1a. You can use code snippets from the labs. Note that the quantitative results are not as important as the qualitative ones. 1. [4 points] Load the dataset and construct the dataloaders for train, validation and test, and visualize the data. Use an 50000 and 10000 images for train and validation, respectively. Which transform is necessary for the samples to be compatible with the models we will create? 2. [6 points] MultiLayer Perceptron (MLP) : Construct a two-hidden-layer MLP with 100 neurons (per layer), ReLU activation functions and a linear output layer to classify FashionMNIST. Train this model for 20 epochs using the cross-entropy loss and the following optimizers : SGD (lr=0.01), SGD with momentum (lr=0.01, momentum=0.9, nesterov=True), Adam (lr=0.01) and Adam (lr=1). Plot the training and validation learning curves (loss against steps and epochs) on a single plot. Comment on your results. 3. [6 points] Convolutional Neural Network (CNN) : Construct a CNN with three convolutional layers (kernel size=3) and 16, 32, and 64 channels, respectively; a non-linearity and one max-pooling layer (kernel size=2) after every convolution; and a final fully connected layer. Train this model with the same four configurations as before, and plot the training and validation learning curves on a single plot. Comment on your results. 4. [6 points] Create a function that computes the number of parameters of a given model. Show the number of parameters for the two models you have used. Does more parameters translate to better performance? Explain. 5. [8 points] PermutedFashionMNIST1 : In this version of the dataset, the pixels are randomly permuted. Visualize the new dataset. Train one MLP and one CNN using SGD with momentum lr=0.01, momentum=0.9, nesterov=True). What do you observe? Is one model more affected than the other? Explain. 6. (Bonus [5 points]) Using any of the tips and tricks seen in class or during the exercise sessions, optimize the validation performance on the original FashionMNIST and report both validation and test performance. Also feel free to explore things beyond the ones seen in class. In any case, explain your decisions. Did you manage to improve the performance? Why?1. Add the tranformation RandomPermutation to the list of transforms. This transform should be the last one.(a) FashionMNIST samples (b) MultiMNIST example (c) Multi-Task model Figure 1 – Samples and models seen in Exercises. 3 Multi-Task Learning with MultiMNIST [30 points] In this exercise you will create a dataset with two tasks from scratch and train a multitask model. Some 1. [7 points] Create a new sample. Implement the function make new sample(x1: torch.Tensor, x2: torch.Tensor) -> torch.Tensor that takes two MNIST samples and creates a new MultiMNIST sample. An example is shown in figure 1b. You can use the following procedure : — Create an empty tensor of size 36 × 36 — Place one image on the top-left corner. Reminder : A MultiMNIST image has dimensions 28 × 28. — Place the other image on the bottom-right corner. For the overlapping pixels take the max value. Visualize five samples with corresponding labels. 2. [5 points] Create the new dataset. Repeat the above procedure 60000 and 10000 times to create the full train and test datasets, respectively. The MultiMNIST train images must be created by samples of the MNIST train dataset. Similarly, for the test dataset. Hint : Use TensorDataset to store the images and their labels. 3. [9 points] Create the Multi-Task model. Adapt the Convolutional Neural Network seen in class to have two outputs by implementing the model in figure 1c. The encoder has two convolutional layers of 10 and 20 channels, respectively, and a Linear layer with 50 output features. Each convolutional layer has kernel size of 5 and is followed by ReLU and Maxpool. Each decoder is a linear layer. 4. [9 points] Train the Multi-Task model. Train the multi-task model and report the accuracy and loss per task. Use the average of task losses as your objective. Hint : Adapt the training and predict functions seen in class to handle two predictions, two losses etc.
This homework involves some coding, with the language of your choice. Present your results with graphs, plots, and data. Submission instruction : The homework must be uploaded on Moodle and the name of the submission must be name surname SCIPER HW2.ipynb. All the results must be included in a single notebook. We recommend you to use Google Colab. 1 Classifying digits with Scikit-learn We are going to use the UCI ML handwritten digits datasets. Using Python and Scikit-learn, they can be downloaded writing : from sklearn.datasets import load digits and then digits = load digits() . This dataset contains ntot = 1797 images, each a set of d = 8 × 8 = 64 pixels, with each pixel being an integer between 0 and 16. The image represents handwritten digits, from 0 to 9. Our goal will be to train a classifier to recognize if the digits are representing even numbers or odd ones1 . 1) Import and prepare the data : — Import the dataset using load digits(return X y=True). — The images are automatically imported as vectors xi,(i = 1,…,n), but so far the labels are numbers between 0 and 9. Change the targets yi to correspond to 0/1 for even and odd numbers respectively. — Split the data into a “training” & a “testing” set (using roughly 60% and 40% ). — Check that these two subsets have roughly the same proportion of labelled numbers. 2) Logistic Regression : — Using linear model.LogisticRegression, train a classifier based on Logistic regression to decide if an image represents an odd or an even digit, comparing `2 and `1 penalty (you can fix the penalty just by setting the corresponding parameter when you call LogisticRegression() to generate the method). — Using GridSearchCV with n splits=5, as seen during the exercise sessions, perform cross validation to fix the optimal hyperparameters (the regularization constant for the penalty) studying the interval λ ∈ [10−1,104] on a logarithmic scale. Then, with the resulting classifier, compute the Accuracy on the Test set, i.e. the rate of correctly classified images. 3) Ridge and Hinge : Repeat these operations for Ridge and Hinge (also called SVM) : — Use linear model.RidgeClassifier to train a classifier based on Ridge, studying λ ∈ [10,107]; — Use svm.LinearSVC to train a classifier based on SVM/Hinge, studying λ ∈ [1,104]. — Compute the accuracy of the classifiers on the Test set as before. 4) Random Forest : Repeat the operations described in point 2 for a random forest classifier : — Use sklearn.ensemble.RandomForestClassifier to train a random forest classifier. Study the tuning of the hyperparameter n estimators with cross validation as in previous points, fix the search interval to be n estimators ∈ [10,104]. — Compute the accuracy of the classifiers on the test set as before. 5) Random Feature : — An alternative is to consider the random feature classifier. This is done as follows : first we transform the vector x containing for each image into another vector u using the following transformation : u )with, σ(x) = 1/(1 + e−x) (1) 1 Here the xi are vectors of dimension d that contain the values of the pixels for a given image, F is a random Gaussian matrix FD×d, where each element is a Gaussian random number of variance 1/d (this can be created using random.normal in numpy). The function σ(·) is applied on each component. This transformation maps all the original vector x ∈ Rd into seemingly ”randomish” vectors u ∈ RD. — Perform the same classification as before with `2-Logistic regression, now using the vector u instead of the vector x. How does the accuracy behave as D increases? (try with D = 2d,D = 4d,D = 8d…). We suggest to use λ ∈ [10−7,10−2] for the CV analysis. 2 Statistical Learning with Nearest-Neighbors 1. First, we assume that we are a “genie” that knows the distribution P(X,Y ). Show that the “best” estimator f(X) in terms of minimizing the expected population risk EX,Y [(f(X) − Y )2] is given by fBayes(X) = E(Y |X). This will serve as our reference point. 2. We are given a set of point x1,…,n and their labels y1,…,n and consider the KNN function (where Nk(X) is the set of k-nearest neighbors of the point X) : (2) i∈Nk(X) Show that the difference between the expected population risk (where we also average over the labels Yi, i = 1,…,n) and the Bayes population risks follows a Bias-Variance decomposition given by RKNN − RBayes Bayes + v (3) with X Bayes(xi) − fBayesand v =(4) k i∈Nk(X) 3. To bound the bias term, we make the following assumptions : (i) we suppose that the regression function fBayes(x) is is L-Lipschitz, and (ii) that the x1,x2,…,xn are are evenly spaced on a d-dimensional unit (hyper-)cube of volume 1. Show that this leads to RKNN − RBayes (5) Pause to think about the relationship between our results and L, k, n, and d. Does this dependence qualitatively make sense? From now on, for the rest of the homework, we shall assume L = σ = 1 for simplicity. 4. To get an idea of what the terms in the bound look like and what the best k might be, plot the individual terms in the bound for the excess risk and the bound itself for n = 100 and d = 1 as a function of k to visualize the trade-off. Find the value of k that minimizes the excess risk. Repeat for different values of n and d. 5. For general n and d, compute analytically the best k, denoted k∗(n,d), as a function of n and d. How does k∗(n,d) changes with d and n? Does this relation with respect to n make sense? 6. Plug this value of k∗(n,d) in the bound (5). Obtain an optimal estimate for the bound on the excess squared risk for k-nearest-neighbor regression in terms of n and d. Pause to think about the result. 7. To understand and get an intuition on the meaning of this result, plot the number of samples n required to achieve excess squared risk of 0.1 as a function of the dimension d. Do you find anything worrying about this plot? 2
This homework involves some coding, with the language of your choice. Present your results with graphs, plots, and data. Jupyter notebooks are a good option, and we recommend you to send your work as a notebook on Google colab. 1 First passage time of a one-dimensional Brownian particle A bit of context : Particles suspended in fluids can be described using the theory of random processes, and in particular Brownian motion. In one dimension, the motion is described by the diffusion equation ∂tp(x,t|x0) = D∂x2p(x,t|x0), where x0 is the particle starting position and D the diffusion constant. It is possible to show, see e.g. Wikipedia, that the (density of) probability that a particle first reaches a point xc (at distance d = |xc − x0| from the starting point) at time t is given by , (1) which is the so-called L´evy distribution. In the following, we will fix the distance d, and starting from the observation of a set of times {ti} we aim to find an estimate of the diffusion constant D. 1. Write explicitly the probability n independent events at times and the normalized log-likelihood . 2. Write a program that simulates n such observations sampled from the probability distribution (1) for fixed d = 2 and D generic. This can be done, for instance, with the scipy.stats.levy method in python, choosing as parameters loc = 0 and scale = d2/(2D). Start with n = 10 observations and plot the (log)likelihood as a function of D. Repeat for n = 20,100, discuss and comment on what you see. 3. We now assume that we are given a set of n observations , without being told the true value D∗. We consider the maximum likelihood estimator ) = argmax (2) D and we shall define the squared error as SE = ( . Create some data sets with n = 10,100,1000 for different values of D∗ ∈ (0,5] and see how the ML estimator performs. Note that finding the maximizer DˆML can be done numerically in python, for instance using the scipy.minimize method. 4. One can show (Bonus 1 : prove it !) that the Fisher information in our problem is given by . (3) Knowing this, another interesting estimator is given by the maximum as posteriori (MAP) with the Jeffreys prior : . (4) Implement this estimator using scipy and repeat the analysis done in point 3. 5. If we average on many realizations (say about a hundred) we can obtain numerically the averaged mean squared error MSE(D∗,D,nˆ ) which is thus a function of n, D∗ and of the estimator Dˆ. Compute and plot, for n = 10,100,1000, the curves MSE(D∗,DˆML,n) and MSE(D∗,DˆJ,n) as a function of D∗. 6. How do the MSE curves at various n compare with the Cram´er-Rao bound for unbiased estimator MSE( (where I(D) is the Fisher information)? How do the Jeffrey and ML estimator behave? Which one would you choose? 7. Bonus 2 : Look at the median of the L´evy distribution on Wikipedia. Can you design an estimator for D from it? Repeat the analysis at point 3 and at points 5 and 6 for this third estimator and compare its performances to the other two. 1 2 Probability bounds and a pooling problem (5) 1. Our starting point is to realize that P(a ≥ b) = P(eλa ≥ eλb) for any λ ≥ 0. Using Markov inequality show that : (6) Hint : Since the Zi are independent, we can write EZ[expPi Zi] = Qi EZi[expZi]. Use the proof strategy discussed in the lectures. 2. Using the value of λ that minimizes the right-hand-side of the former equation, show that with Hint : Instead of differentiating the right-hand-side of (6), one can differentiate its log, since the logarithm is monotone and thus its extrema are the same of its argument. 3. Show that , and that 4 for any s.t. 0 . 4. Use Taylor’s theorem (that states that 2 for some unknown ˜, and where the prime stands for derivative with respect to ) to show that , and prove the inequality (5). Similarly, it is possible to show (Bonus 3 : prove it) that so that (7) 5. The most important use of such a bound is in terms of pooling problems. Suppose you want to know what fraction of the population in a country approves its current president : how many people should you ask to be confident, with probability at least 95 percent, that the error in estimating the fraction of people who approves the president is correct within one percent (so that ˆp is in [p − 0.01,p + 0.01] with 95% probability)? 6. Compare the number m∗ you find this way with what you observe when performing numerical experiments in python : — Define a function that takes the number of people m and the probability p as arguments and returns a random array of m votes. Hint : You can generate Bernoulli-distributed samples in python using the scipy.stats.bernoulli method. — Starting with fixed values of m = m∗ and p ∈ {0.2,0.5,0.8} — Use this function evaluated in m∗ and p to simulate polls. — Just by averaging the generated votes, estimate p. — Quantify the probability that ˆp is correct within one percent. — Which values of p seem to be harder to estimate? Do you find that the bound is accurate, or does it grossly overestimate the needed number? — For each p, repeat for different values of m to find the value that (more or less) gives an estimate which is correct within one percent with 95% probability. — Bonus 4 : Plot the behaviour of the probability of error as a function of p for values of m ∈ [10,104], and compare it with the theoretical m∗. 2
Display InterfaceLab Objectives In this lab, we will learn how to display images on a VGA terminal using the DE1-SoC Computer’s videoout port and implement a line-drawing algorithm. Introduction The DE1-SoC Computer includes a video-out port with a VGA controller that can be connected to a standard VGA monitor. The VGA controller supports a screen resolution of 640 × 480. The image that is displayed by the VGA controller is derived from two sources: a pixel buffer and a character buffer, however, we will only discuss and use the pixel buffer in this lab. Pixel Buffer red green blue 0…0 y x 0 The pixel buffer for the video-out port holds the data (color) for each pixel that is displayed by the VGA controller, with the coordinate ( , ) referring to the top-left corner of the image (Figure 1). 15 ⋯ 11 10 ⋯ 5 4 ⋯ 0 (a) Pixel color 31 ⋯ 18 17 ⋯ 10 9 ⋯ 1 0 (b) Pixel (x,y) offsetFigure 2: Data layout of pixel color values and address offsets. The color of a pixel is represented within a 16-bit halfword, with five bits for the red and blue components and six bits for the green component (Figure 2a). The address of an individual pixel in the pixel buffer is given by the sum of a base address and an (,) offset (represented as seen in Figure 2b). In the DE1-SoC, the default base address of the pixel buffer is , which corresponds to the starting address of the FPGA on-chip memory. • e.g., pixel ( , ) has address pixel ( , ) has address pixel ( ,1) has address pixel ( , ) addr: You can create an image by writing color values into the pixel addresses as described above. A dedicated pixel buffer controller reads this pixel data from the memory and sends it to the VGA display. The controller reads the pixel data in sequential address order, starting with the pixel data corresponding to the upper-left corner of the VGA screen and proceeding horizontally, then vertically until it reaches the data for the lower-right corner; this process repeats continuously. You can modify the pixel data at any time by writing to the pixel addresses. Writes to the pixel buffer are automatically interleaved in the hardware with the read operations that are performed by the pixel buffer controller. Drawing Drawing a line on a screen requires coloring pixels between two points (x1,y1) and (x2,y2), such that the pixels represent the desired line as closely as possible. For example, if we want to draw a line between points ( , ) and ( , ), we cannot draw the line precisely – we can only draw a shape that is similar to the line by coloring the pixels that fall closest to the line’s ideal location on the screen. This is illustrated in Figure 3, where the squares represent pixels on the screen.We can use algebra to determine which pixels to color based on the slope and end points of the line. The slope of our example line is (2 − 1)/(2 − 1) = 4/11. Starting at one endpoint , we move along the -axis to the other endpoint and compute the -position for each -position as: = 1 + slope × ( − 1). Thus, for = 2, which you can think of as a column of pixels, the -position of the pixel for the line is 1 + (4/11) × (2 − 1) = 15/11. Since pixel locations are defined by integer values, we round this result to the nearest integer, and determine that we should color the pixel at . For = 3, we compute = 19/11 and color the pixel at . Similarly, we compute pixels to color through = 12 = 2. This approach (compute a -position for every -position) has drawbacks when a line is steep (i.e., |slope| > 1), since it will span more rows than it does columns. In this case, our calculations will not produce a smooth-looking line and will not work at all for a vertical line. To address this problem, we can alter the algorithm to instead move along the -axis when a line is steep. With this change, we can implement what is known as Bresenham’s line algorithm (two different pseudocodes are given in Figure 4 and the first link on the next page – either will work for this lab). In Figure 4, the first 15 lines of the algorithm make the needed adjustments depending on whether or not the line is steep. Lines 17 to 22 then increment the variable one step at a time and computes the corresponding -value. The -value is incremented when needed to stay as close to the ideal location of the line as possible. Bresenham’s line algorithm calculates an variable to decide whether or not to increment each -value. The version of the algorithm shown in Figure 4 uses only integers. This arithmetic is signed, so pay attention to your data types and operators in SystemVerilog and test both positive and negative values, where appropriate. If you would like to learn more about how Bresenham’s line algorithm works: • http://members.chello.at/easyfilter/bresenham.html • https://en.wikipedia.org/wiki/Bresenham%27s_line_algorithm 1 draw_line(x0, x1, y0, y1) 2 3 boolean is_steep = abs(y1-y0) > abs(x1-x0) 4 if is_steep then 5 swap(x0, y0) 6 swap(x1, y1) 7 if x0 > x1 then 8 swap(x0, x1) 9 swap(y0, y1) 10 11 int deltax = x1-x0 12 int deltay = abs(y1-y0) 13 int error = -(deltax/2) 14 int y = y0 15 if y0 < y1 then y_step = 1 else y_step = -1 16 17 for x from x0 to x1 18 if is_steep then draw_pixel(y,x) else draw_pixel(x,y) 19 error = error + deltay 20 if error >= 0 then 21 y = y + y_step 22 error = error − deltaxFigure 4: Pseudocode for Bresenham’s line algorithm. Task #1 – VGA Connection (no deliverables) 1) Download from Canvas (Files → Labs → lab5). 2) Unzip the Quartus project folder, which contains some files that will help you work with the VGA output of the DE1-SoC. In particular, there are three SystemVerilog files: a) is a driver for the VGA port of the board. You don’t need to edit or understand this file, but you might notice that it uses a 38,400-byte framebuffer register, similar to what was described above. The ternary operator on the last line of this file controls the colors of the lines you’ll be drawing (only black or white for this lab). b) is a skeleton file for you to implement Bresenham’s line algorithm. c) is a top-level module that instantiates both of the above modules. This should compile and function as-is, but you will modify it for later tasks. 3) Verify that the project compiles on LabsLand and produces an output on your monitor. a) On LabsLand, choose the “VGA” option for the interface (shown in Figure 5). b) Once the code is synthesized and uploaded onto an FPGA, the VGA monitor will initially produce a rainbow gradient before turning black (shown in Figure 6).Task #2 – Implement the Line-Drawing Algorithm Some notes about the • It takes inputs , corresponding to the coordinate pairs • On positive edges of the input clock , the outputs and form a coordinate pair on the line between and . On any given clock cycle, and should each change by at most one pixel, since we’re guaranteed to use |slope| ≤ 1. • As indicated in the file, you’ll need to create some local registers (i.e., variables) to keep track of things. Notice that the example signal is declared as and is a bit longer than the coordinate ports. Your task is to draw a line between any two arbitrary points on the monitor. However, Bresenham’s line algorithm can get complicated and you need to handle drawing to the left, right, up, or down, and both steep and gradual lines. Instead of doing this all at once, you’ll probably want to work in smaller steps. The following are suggestions on how to approach this problem, but you can complete this task in whatever way makes the most sense to you. 1) Set = (in case of monitor cut-off) and use to draw horizontal lines – i.e., implement the for-loop mechanics and > conditional. 2) Modify your algorithm to handle perfectly diagonal lines from any arbitrary starting point – i.e., handle 3) Modify your algorithm to handle lines with gradual slopes, e.g., from – i.e., handle 4) Modify your algorithm to handle lines with steep slopes, e.g., from – i.e., handle While you won’t be demoing this task, you will need to include simulation results for drawing left-up, right-up, left-down, right-down for both steep and gradual slopes. Task #3 – Animate an Object Modify to implement the following: 1) Use your line algorithm to draw a line on the monitor and animate it to move around the screen. a) You can be as creative as you want with this animation, but make sure that it demonstrates the full functionality of your line-drawing algorithm. You should draw lines to the left and right with negative, positive, steep and shallow slopes as well as horizontal and vertical lines. 2) Implement a reset that, when activated, clears the screen by drawing every pixel to be black. a) You’ll need to modify the port connection to the module to choose between drawing black or white. Note that you don’t need to provide simulation results for your animation, but instead will include a drawing or image of what it produces.Lab Demonstration/Turn-In Requirements In-Person Demo • Demonstrate your working Task #3 that animates an object moving around the screen. • Demonstrate your working Task #3 reset that clears the monitor. Lab Report (submit as PDF on Gradescope) • Include the required Design Procedure, Results, and Experience Report sections. • Don’t forget to also submit your SystemVerilog files ( ), including testbenches!Lab 5 Rubric Grading Criteria Points Name, student ID, lab number 2 pts Design Procedure ▪ A diagram must be included for your algorithm. ▪ The design of your animation algorithm must be described and please include a diagram, if appropriate. 14 pts Results ▪ Correct behavior must be shown for drawing left-up, right-up, left-down, rightdown for both steep and gradual slopes. ▪ Include a drawing or image of your final animation. 18 pts Experience Report 6 pts SystemVerilog code uploaded 5 pts Code Style 5 pts LAB DEMO 80 pts ▪ Bonus points available for particularly impressive animations. (10 pts) 130 pts
